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SUMMARY 


This report volume presents a unified theory for aerodynamics and noise of 
single rotation propellers. Aerodynamic topics include calculation of perfor- 
mance, blade load distribution, non-uniform wake flow fields, and non-uniform 
upstream (potential) flow fields. Blade loading can be steady or unsteady due to 
fixed distortion, counter-rotating wakes, or blade vibration. The aerodynamic 
theory is based on the pressure potential method and is therefore basically lin- 
ear. However, non-linear effects associated with finite axial induction and blade 
vortex flow are included via approximate methods. Acoustic topics include radia- 
tion of noise caused by blade thickness, steady loading (including vortex lift), 
and unsteady loading. The loading orientation includes the radial component. 

The derivation begins with the development, via the acoustic analogy, of a 
general expression for disturbance pressure as an integral over the blade surface 
of the loading and thickness distribution. This immediately provides expressions 
for near field noise predictions. These expressions are specialized for far field 
noise calculation using the method of stationary phase. A section on sound power 
and wave drag is presented to evaluate the contribution of acoustic wave energy to 
aerodynamic performance loss. For unswept, supersonic tip speed blades, the 
acoustic power cannot be neglected in comparison with the shaft power. 

For steady aerodynamic use, the integral formula for pressure is converted to 
a formula for downwash at control points on the blade surface as a function of 
loading and thickness distributions. Since the downwash is considered known from 
the blade camber and angle of attack distributions, the integral equation must be 
inverted to find the blade loading. This is accomplished in a manner similar to 
that of wing methods by discretizing the loading, converting the integral equation 
to a matrix equation, and inverting the matrix. Induced drag (or vortex drag) is 
computed in the Trefftz plane. A section is devoted to showing the relationship 
between this new lifting surface theory and classical propeller momentum theory. 
When the blade loading has been calculated, the wake and upstream potential field 
calculations are done using velocity formulas similar to the downwash equations 
but for field points off the blade. 

For unsteady aerodynamic application, the downwash is considered known from 
the blade vibration. The integral equation is solved for unsteady loading by the 
same inversion scheme described above. The only numerical results presented in 
this volume are some verification cases for the 3D unsteady loading theory. 
Calculations of performance, wakes, and noise are presented in Volume III. 



SECTION 1 
INTRODUCTION 


This report presents results of program for the generation of an accurate, effi- 
cient computer prediction code for noise of advanced, single rotation, turbo- 
props (Prop-Fans) such as the SR3 model shown in Figure 1.1. The code is based 
on a linearized theory developed at Hamilton Standard in which aerodynamics and 
acoustics are treated as a unified process. Both steady and unsteady blade 
loading are treated. Capabilities include prediction of steady airload distri- 
butions and associated aerodynamic performance, unsteady blade pressure response 
to gust interaction or blade vibration, noise fields associated with thickness 
and steady and unsteady loading, and wake velocity fields associated with steady 
loading. The code was developed on the Hamilton Standard IBM computer and has 
now been installed on the Cray XMP at NASA-Lewis. 

The work had its genesis in the frequency domain acoustic theory developed at 
Hamilton Standard in the late 1970' s. It was found that the method used for 
near field noise predictions could be adapted as a lifting surface theory for 
aerodynamic work via the pressure potential (or acceleration potential) tech- 
nique that has been used for both wings and ducted turbomachinery. In the first 
realization of the theory for propellers (prior to the contract) , the blade 
loading was represented in a quasi-vortex lattice form. Under the contract, 
this was upgraded to true lifting surface loading. Originally, it was believed 
that a purely linear approach for both aerodynamics and noise would be adequate. 
However, in the course of the contract, two sources of non-linearity in the 
steady aerodynamics became apparent and were found to be a significant factor at 
takeoff conditions. The first is related to the fact the steady axial induced 
velocity may be of the same order of magnitude as the flight speed and the sec- 
ond is the formation of leading edge vortices which increase lift and redistrib- 
ute loading. The contract was amended to deal with both of these phenomena. 

The final report is divided into 5 volumes as outlined in the Abstract. This 
volume (Volume I) gives the theoretical derivations and background for the aero- 
dynamics and acoustics of rotating, unducted blades. Verification by comparison 
with test data is the subject of Volume III, however, this volume does include 
some verification by comparison with related theories. In particular, the rela- 
tionship to classical propeller theories is established and treatment of high 
order singularities is compared with that in an analogous compressor method. 

The volume is divided into 13 sections, each depending on previous material. 
However, for readers Interested is specific topics, not all of the sections need 
to be understood. A guide to the report sections is provided on the next page 
to help in understanding the flow of the derivations. 
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SECTION 2 
BACKGROUND 


This section presents a brief overview of some of the analytical methods that 
can be used for propeller aerodynamics and acoustics. Advantages of the pressure 
potential method are explained and some of its relevant history is reviewed. 

There is a broad spectrum of methods that could be applied to the propeller 
problem and each has advantages ranging from economy on one end of the scale to 
accuracy at the other end. At the top end of the scale, the Euler and Navier- Stokes 
methods that are currently under development for Prop-Fans may, at some time, be 
used routinely for propeller acoustics. Indeed, Euler aerodynamic methods are 
currently being used to provide airloading as input for acoustic theories (ref. 1) 
with considerable success. The main drawback to this approach is the long computer 
running time for fully converged aerodynamic results. Calculations for steady 
airloading are rapidly becoming manageable but it appears that affordable unsteady 
propeller aerodynamic calculations with these fully numerical schemes are still some 
time in the future . 

The use of an Euler calculation for a direct acoustic prediction might seem 
like a tractable approach for the near field. However, attempts so far (ref. 2) 
have suffered from poor grid resolution. Grid points must be packed closely near 
the blade surface for good loading results. If shock waves are to be followed out 
into the acoustic field, then close packing is needed there also. The result is a 
requirement for more grid points than are currently possible. 

Potential field methods are based on less complete versions of the fluid 
dynamic equations but are of interest because powerful solution techniques are 
available. In particular, Green's function methods provide solutions for all space 
without the requirement for finite volume meshes. They deal with radiating waves 
naturally and can thus be used successfully for both aerodynamics and acoustics. 
Potential equations fall into the categories of velocity potential and pressure (or 
acceleration) potential. The major difficulty with velocity potential theory is 
that sources must be placed on both the blades and on their wakes. This leads to 
extra computational complexity, particularly for Prop-Fan cruise conditions with 
subsonic roots and supersonic tips. With the pressure potential formulation, how- 
ever, sources need only be placed on the lifting surfaces and not on the wakes. 
Furthermore, the mixed subsonic, supersonic condition can be handled with Fourier 
transform methods. The pressure potential does lead to higher order singularities, 
but these can be handled analytically resulting in a practical computational method. 

Pressure potential methods have been in use for many years for wing steady and 
unsteady aerodynamics and are treated in standard text books. For example, Bisp- 
linghoff , Ashley, and Halfman (ref. 3) explain the method in some detail for both 
incompressible and compressible flow. Also, Bisplinghoff and Ashley (ref. 4) reduced 
the method to operator form for aeroelastic applications. In simplest terms, the 
pressure potential method works as follows. An analytical expression is found for 
the pressure field p around the wing or body of interest. The gradient of p is 
then used in the momentum equation 


st - h? - v p> 
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where the left side represents the acceleration of a fluid particle and F is the 
body force, if any, acting on it. To find velocity, the acceleration of the par- 
ticle is integrated from some place where the velocity is presumed known, typically 
upstream infinity. Body forces are avoided by dealing with particles that don't 
quite come into contact with the wing. (Flow at the wing is handled by a limiting 
process.) Linearization enters the problem in 3 significant areas. First, the 
pressure field is usually found by a solution of the linear wave equation. This 
means that bow and trailing shocks are approximated by Mach waves; attached, mid- 
chord shocks are not represented. The second linearization is in setting the den- 
sity p to its ambient value p°. This is justified because the above equation is 

order in the disturbance pressure so that disturbances in density would bring 
in second order quantities. The third entry of linearization is in approximating 
the particle trajectory by its undisturbed path. This seems to be satisfactory for 
wings because there is little fluid deflection until the particles are very near the 
leading edge. (The same is not necessarily true of propellers.) The pressure 
potential method has generated countless papers and computer programs for wing 
computation. The paper by Landahl and Stark (ref. 5) gives an excellent review of 
the subject. 

The pressure potential method has also been applied to ducted rotors. The 
analyses by Namba (ref. 6) and by Lordi and Homicz (ref. 7) follow the general 
recipe given above, but with vast differences in the details. Namba' s theory was 
used for aerodynamics and acoustics, whereas Lordi and Homicz applied theirs only to 
computation of blade loading. Their work is particularly interesting, however, in 
that it brings new insight on the role of body forces and high order singularities. 
These analyses use Green's functions for constant diameter, hard walled ducts and 
thus do not permit the kind of contracting inflow that characterizes propellers . 
Instead they do produce a pressure jump at the rotor that persists to downstream 
infinity. 

For propellers, the earliest application of the pressure potential method seems 
to be that of Rondo (ref. 8) who worked out a version of the steady flow equations 
during the World War II. His expressions were not written out explicitly and appar- 
ently have never been coded for computers. Much later Tsakonas and his co-workers 
at the Stevens Institute of Technology developed an incompressible pressure poten- 
tial method (ref. 9) for marine propellers that has been used for steady performance 
and for unsteady blade loading. Their integral equation is equivalent to that of 
the present report in the limit of zero Mach number. 

The evolution of the present work has been as follows. A helicoidal surface 
theory (ref. 10) was derived for prediction of propeller noise caused by thickness 
and steady loading. This was written directly in the time domain and worked quite 
well for ordinary propellers even with supersonic tips. However, the combination of 
swept blades and supersonic tips lead to waveform predictions with a jagged charac- 
ter caused by numerical difficulties associated with retarded source location and 
numerical differentiation. These numerical problems were overcome by changing to a 
frequency domain approach which eliminated the need to deal with retarded source 
calculations. The first 2 papers using this method were restricted to the far field 
and to steady loading (refs 11,12). Near field formulas were then worked out, again 
for steady loading (refs. 13) . The far field unsteady loading theory (ref. 14) 
next was presented to deal with unsteadiness caused by interaction between rotors of 
a counter-rotation propeller. Also, near field equations for noise caused by 
unsteady loading were derived but not coded or published. 

At the same time, it was discovered that the near field equations were suitable 
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for use in pressure potential aerodynamic analysis. The resulting theory (ref. 15) 
was presented showing the connection between the aerodynamics and acoustics and how 
the new theory related to previous theories for propeller noise and hydrodynamics 
and to theories for wing steady flow and flutter analysis. The first application of 
the blade loading theory was given in Reference 16. This treated the steady loading 
problem only and represented the loading in a sort of vortex lattice scheme. 

That was the status of the analysis when the current contract work began. The 
major extensions of the work planned under the original contract included upgrading 
the loading representation from vortex lattice to true lifting surface, coding the 
unsteady loading noise and airload theories, adding wake prediction capability, and 
bringing all of the code together in a unified computer system. The work applies 
strictly to single rotation propellers although the wake and unsteady lift response 
capability needed for counter-rotation are included. Also, only linear theory was 
to be included. This was reasonably satisfactory for the high speed cruise cases 
where axial induction is small. However, for takeoff cases, it was realized that 
the large axial induction (compared to the flight speed), violated the linearization 
assumption described above that permitted integrating the fluid particle acceler- 
ation along its undisturbed path. This problem does not occur in textbook wing 
problems nor in ducted rotor problems where contraction ahead of the rotor cannot 
occur. It was not anticipated in the propeller problem. The contract was amended 
to address this problem and a relatively simple scheme was developed. In the pro- 
cess, the connection between the new theory and classical propeller momentum theory 
was established so that the continuum of methods can be now be understood. 

The other non-linear effect that was not anticipated was the appearance of 
leading edge vortices at takeoff conditions. These correspond to the well known 
vortices on delta wing aircraft that are helpful in providing extra lift for takeoff 
and landing at low speeds. The conditions necessary for their formation are sharp 
leading edges, sweep, and high loading, all of which are present on Prop -Fans at 
takeoff. It turns out that the presence of leading edge vortices has been known on 
marine propellers for some time. Kerwin's review article (ref. 17) includes a 
photograph of a leading edge vortex made visible by cavitation bubbles. A sketch 
showing their appearance in a water tunnel experiment was presented in Reference 18 
and is reproduced In Figure 2. The similarity to typical Prop-Fan geometry is 
remarkable. It seemed very probable that the same phenomenon was occurring on 
Prop-Fan blades and some flow visualization experiments were initiated in the United 
Technologies Acoustic Research Tunnel to demonstrate their presence. The first 
successful results were via surface oil flow patterns made visible with florescent 
dye as shown in Figure 3. The streaks show the oil flow to be predominantly in the 
radial direction under the influence of centrifugal force. However, there Is a 
region were they change direction, as shown in the inset to the figure. This is the 
region of re -attachment of the leading edge vortex. The wind tunnel experiment was 
expanded to document the flow patterns over a wide range of operating conditions and 
the results were released in 1987 (ref. 19). Similar results have now been produced 
in experiments at NASA-Lewis (ref. 20) and there now seems to be little doubt about 
the presence of the vortices . 

Since the leading edge vortices tend to produce extra lift and redistribute the 
air loading, the contract was also amended to deal with them. The treatment is via 
the suction analogy of Polhamus (ref. 21) and Lamar (ref. 22) in which leading edge 
thrust, which can be computed with a lifting surface panel method, Is used to pre- 
dict the extra leading edge loading associated with the separation and reattachment 
of the vortices. Application of the suction analogy to propellers is described in 
Volume III of this report. 
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SECTION 3 

GENERAL PRESSURE THEORY 

This section presents a derivation of the general equation for the pressure 
field of a propeller. The theory is based on the linear wave equation with heli- 
cally convected thickness and loading sources representing the action of the blades 
on the air. A solution is achieved using the free -space Green's function and 
Fourier transforms . The general equation is needed in later sections where working 
equations that have been coded for near and far field noise prediction are devel- 
oped. It is also the starting point for the performance analysis in Section 7. 


Differential Equation 

The linear governing equation for pressure disturbances is 

v2 P - \ —2 - -S(r.t) (3-D 

c o at 2 

where the source function is given by 

S(r.t) = p - V.f (3.2) 

— o at ~ 

Equation 3 . 1 has solution 

rr £ (t-r-R/c ) 

P(r,t) = JJsCr^r) ^ dr^dr (3.3) 

where the distance between the field point r and the source point r 0 is 

R - |r - r Q | (3.4) 

The equations above can be found in many standard references including, for 
example, Goldstein (ref. 23). In the source term, q is the volume displacement/unit 
time/unit volume and f is the force/unit volume acting on the fluid. In Gold- 
stein's book q represents an actual volume injection but here we use it to repre- 
sent the volume displacement effect of the blades; the resulting pressure waves are 
traditionally called thickness noise. Similarly, f was derived as a volume source 
but is used herein to represent the blade loading in terms of surface pressure. The 
associated radiation is called loading noise. 

In Equation 3.3, 5/4 jtR is the free space Green's function representing spher- 
ically propagating wavelets from the source . The double integral sums the source 
contributions from all space and time and the integration limits, which are not 
specified in Equation 3.3, cover all contributing regions. The problem to be 
treated here involves blades that are translating and rotating through the fluid. 
Other treatments of related problems use coordinate systems translating with the 
propeller (the aircraft system) or both translating and rotating with it (the blade 
system). Equation 3.1, however, is written for a coordinate system fixed in the 
fluid so that, at the beginning of the analysis, we will be dealing with the tran- 
sient problem of a propeller fly-by. Because of the Fourier methods used, it will 
be seen that this leads to a simplified theoretical development and, in the final 
pressure equation, the field point may be moved to either the aircraft coordinate 
system or the blade coordinate system via a trivial transformation. Furthermore, 
using this approach, the only Mach number that appears is that associated with 
forward flight, which applies at all stations on the blade. Blade fixed systems 
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have to deal with varying Mach number along the radius and a change in the type of 
differential equation for supersonic section speeds. This requires some sort of a 
patch in solutions at the sonic radius. 

Previous treatments by the author of the propeller aero-acoustic problem were 
based on the acoustic analogy equation, which can also be found in Goldstein (ref. 
23). This brings in the surfaces explicitly; however, it was found when the bound- 
ary conditions are linearized through the thin blade approximation, that the dis- 
tinction between the acoustic analogy and Equation 3.1 vanishes. Thus, the simpler 
approach is taken here. 


Source Representation 

Before Equation 3.3 can be integrated, the source term S will be represented 
in propeller terminology starting with the sketch in Figure 4. On the left is shown 
a cylindrical coordinate system centered on the propeller axis of rotation with the 
pitch change axis (PCA) coincident with x=0 , <f>~0 at time t=0. </> is positive in 

the direction of rotation regardless of the direction of rotation. The sketch at 
the top shows a blade (say blade number 0) intersected by a cylindrical surface of 
radius r. If this surface were cut and rolled out flat, it would have the appear- 
ance of the sketch at the bottom in Figure 4. It is most convenient to formulate 
the source description in coordinates 7, £ aligned with the direction of advance of 
the airfoil section at each radius. Thus, the helicoidal surface (£=0) swept out 
by the pitch change axis advancing at speed V and rotating at angular speed Q is 
analogous to the reference plane used in wing lifting surface methods. Once the 
source is formulated in these coordinates, the description will be converted to 
conventional cylindrical coordinates before the integration of Equation 3.3. This 
will avoid any problems which might arise from the use of the 7 ,|,r system, which 
is non- orthogonal . 


Thickness Source - At time t=0, the rate of volume displacement per unit 
span is given by Uh' ( 7 ) . h Is the airfoil section thickness distribution shown in 
Figure 4 at radius r and U Is the section helical advance speed. The prime on h 
denotes differentiation with respect to the argument. At times other than t=0, 
motion of the source in the negative 7 direction results in an argument shift 
yielding Uh'( 7+Ut). The line of action of the source Is placed on the offset 
surface via the delta function 6(£+FA) where FA is the face alignment. Thus, the 
rate of volume displacement can be written 

q = 6(£ + FA) Uh' ( 7+Ut) (3,5) 

or equivalently 

q - 8(Z + FA) ± h( 7 +Ut) (3.6) 

so that the complete thickness source description in 7,( coordinates is 

P *W~ P , S<( + FA) ^3 h <T fUt > (3-?) 

The thickness distribution Is represented in terms of its transform via the follow- 
ing Fourier transform pair. 
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(3.8) 


fja 

V> v (k 7 ) - f h( 7 )e 7 d 7 


M7) - w «' U, X < 3 - 9) 

Inserting Equation 3 . 9 into Equation 3 . 7 and performing the time derivative gives 

3q 


3t 


P U . -ik (i+Ut) 

-fa- 5(^+ FA)/k 2 ^(k 7 )e 7 dk 7 


(3.10) 


To express Equation 3.10 in cylindrical coordinates, consider the r— constant surface 
to be unwrapped onto a plane. Then the 7 , £ system can be rotated into x,$ using 
the relations 


V 

’U 


x 


QE r0 

U 9 


> Qr V j 

^ = IT u ^ 


where the coefficients V/U and fir/U are tied to the advance triangle shown in 
Figure 4. The inverses of these equations are shown for future reference 


** - - if T - U* 


V . Qr c 

x - - u 7 + tt c 


(3.11) 

(3.12) 

(3.11a) 

(3.12a) 


For the delta function in Equation 3.10, we first note that it is even with respect 
to argument and then substitute Equation 3.12 for f 


*(£+ FA) - «[¥£(* - ^ - feFA)] 


Vr 


^ - T ■ W FA) 


(3.13) 


The second line in Equation 3.13 follows from the integral properties of any func- 
tion, say f (ax) : Jf(ax)dx = 1/a Jf(x)dx. Now, with the substitution of Equation 
3,11 for 7 , the thickness source has the form 


5q ' P ° V At A ft* U FA x 

- = WE 5( ^ ' T ' w FA) 


at 

x Jk^ v (k 7 ) e 1 ’ 7 ' 8 *' u 


.. ,V nr .. 
ik (-x+-rrr<£> 


-ik Ut 
e 7 dk 


(3.14) 


But the source exists only on <^=Ox/V + (U/Vr)FA so that this can be substituted 
for <f> in the exponential. The result, after recognizing that 


defining 


and 


U 2 - V 2 + n 2 r 2 

(3.15) 

<7= U/V 

(3.16) 


(3.17) 

the delta function via 


6(a) - 1 e in “ 

Z7r n= -® 

(3.18) 


- 9 - 



is 


P 


o 


3q 

dt 


4?r 2 r 


i {*** 

n=-« % 




. nfl v 

lCK x” "v~ )x 


i _ nU NTA 

a ( U K x Vr )FA 


-iK Vt 


dKx} 


(3.19) 


where K* is the ordinary x wavenumber. Equation 3.19 is the complete space -time 
description of the rotating thickness source in cylindrical coordinates. It can be 
verified that Equation 3.19 is time independent in the blade coordinate system 
x i> < t> l by substituting x=x x +Vt and ^-^4-Ot and verifying that time drops out of 
the right hand side. 


Loading Source - The effect of loading is handled in similar fashion. Force 
per unit volume f. is expressed in terms of the surface force per unit area 

fj according to 


f. - - 6 (£+FA) f ,( 7 +Ut)e' i ‘'° t (3.20) 

where, again, the delta function places the load on the blade surface, j- 1 , 2,3 
denotes the x,^,r force components and the exponential accounts for the fact that 
the force may be harmonic at frequency The significance of the minus sign in 

Equation 3.20 is that f ^ represents force on the fluid and fl represents force on 

the blades. The force transforms are defined in parallel with the thickness trans- 
form in Equations 3.8 and 3.9. 


V-.(V “ ft (7) 


ik 7 

e 7 d 7 

(3.21) 

-ik t 

e 7 dk 7 

(3.22) 


and the same operations applied to the thickness source lead to 


x e 


4?r 2 r 

i ({ jfVvr )FA ' iK - vt 




By applying the divergence operator in cylindrical coordinates, 


V.f 


^2 

3x 


+ F 


d<t> 


dr 


(rfj 


the loading source becomes 



, ,nr„ ™1 \ . nfi . 

^ ¥«.*>» J b(K x ~)X e - i(Kx V + Vo,t dKx 




(3.23) 


(3.24) 


(3.25) 


The source function is put into final form for integration by shifting the 
origin of the variable by w Q /V and defining the unsteady loading order q according 
to 


H> = qft (3.26) 

The result, including thickness and loading effects, expressed in source coordinates 
r o> *o’ x o> r > is 
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where 


S(r 0 ,r) = 


in# 


47T 2 r „ 


■j>„( 


K x -q(l/V 'j 

% J 


i(K x - 


n+q 


)x„ 




(3.27) 


and 



(3.28) 


(3.29) 


is the phase shift due to face alignment or offset of the blade section from the 
pitch change axis . 


Integration of Source Function 

We are now ready to insert the source function, just derived into the Green's 
function integral, Equation 3.3, and to perform the r, x 0 , and <f> Q integrations noting 

that the volume element is 

dr = r d<j> drdx (3.30) 

The t integral can be done immediately. The result, after shifting the origin of 
the x integration by x, is 


P(l>t) 


16 tt 3JJZj 


_ n+q 


where 




J Jr^+x? 

-eo ^00 


in which we have defined 


with 


r = K + *o 


= r 2 + r 2 - 2rr 0 cos 

Integral I is evaluated in the appendix to this section with the result that 


V i ^ V ^ ^ o 


i - i* l e 


(3.31) 

(3.32) 

(3.33) 

(3.34) 

(3.35) 


where J is an ordinary Bessel function, is a Hankel function of the first 

kind, and K r is the radial wavenumber. The notations r < and r > denote the 
lesser or greater of r and r 0 . 




(3.36) 
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(3.37) 


with branch cuts shown in Figure 5 . K r can also be expressed as 

K r = | K r | e i0 


where the angle 


6*+6 ±n 


2 


(3.38) 


and the branch points are 


k + = 


(n+q)£VV 


and 


(n+q)fl/V 


1+M* ““ “ 1-f^ 

It will be shown later that these correspond to outgoing waves . 


(3.39) 


We now return to Equation 3.35 and substitute it into Equation 3.31. The <f> 


integral is found to be 


f * 1 


(n-m)$ 


° d 4 > - 2ttS 
T o ran 


(3.40) 


The Kroneker delta enables the m- summation so that the pressure simplifies to 

p(r.t) - jl I e« rr e- ,K < Vt 

on n flco r h-« ^ ° J 


x J„(r«K r )H^ 1> (r>K r ) dl^dr,, (3.41) 

In Equation 3.41, the arguments of the J and H functions are complex. For 
numerical evaluation it is more convenient to deal with real arguments. With the 
branches shown in Figure 3 the JH combination can be interpreted as follows for K* 
on the real axis. 


2 9 

Whenever K r < 0, including the case where n+q=0, then 

J n( r < K r)^ n ( r > K r> “ W t ( r < I K r I * ^ (r> I * r I ) (3.42) 

where I n and are modified Bessel functions. In this range of K^, the wave 
pattern decays exponentially with increasing radius. Note that K's with various 
subscripts are used here with different meanings. To assist in interpretation, they 
are defined when they first appear in the text and can also be found In the symbol 
list given in Section 16. 

For > 0, waves are propagating and we use 

J n (r < K s )H^ 1, (r ; ^ r ) - J n ( r< |K r | )H^ 1) (r > |K r | ) for n+q>0 (3.43) 

and 

J n( r < K r)^ n ( r > K r) “ - \ ( r < I K r I ) H^ Z) (r, | K r | ) for n+q<0 (3.44) 

This Is as far as the integrations can be carried analytically. The remaining 
wavenumber integration, radial integration, and harmonic summation must be performed 
numerically. As will be seen later, the shapes of the chordwise pressure and thick- 
ness distributions will be chosen so that the source transform v^ n can be pre- 
computed and stored in numerical or analytical form. 

Equation 3.41 is the pressure field of a one-bladed rotor represented in the 
r,$,x coordinate system fixed in the fluid. Before proceeding to deal with multi- 
blade effects, we will demonstrate the change in form for observers in the aircraft- 
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fixed frame and in the blade -fixed frame. The x, <j> , t information is contained 
exclusively in the exponential combination 

E - n0 + (k, - ^ 0]x - K*Vt (3.45) 

Most of the acoustic analysis that follows will be for an observer in the 
aircraft - fixed coordinate system. His axial position x^ with respect to the plane 
of rotation is related to x by 

x=x 1 +Vt (3.46) 

With this substitution, the exponential becomes 

E = u<f> - (n + q)Clt + K* - ^ «]x x (3.47) 

Here, the transient effect, I^Vt has vanished so that the translating observer 
senses only periodic components with frequencies given by (n+q)ft, as required from 
physical considerations. If the observer is now moved onto the rotating blade via 
the transformation 

0-^+Ot (3.48) 

then the exponential becomes 

E - nfj - qOt + [k* - ^n)x 1 (3.49) 

so that the only frequency sensed by the blade- fixed observer is the source 
frequency qft. 


Normalizations 

Up to this point, the development has been in terms of dimensional variables 
because they provide a better physical understanding of the problem. However, 
computer program code and input are best written in non-dimensional (or scale- 
independent) variables. Conversions of the radiation equation and source functions 
to non-dimensional terms are handled in this section. 


Source and field point radius ratios are defined as 

z o=r 0 / r T and z=r/r T (3.50) 

and ratios of tip speed and section helical speed to flight speed are given by 

a-nr T /V, cr=U/V=*Jl+a z z 2 , a 0 = u 0 (3.51) 

Also, the wavenumber integration variable is rescaled using 

K x r I =au (3.52) 

With these changes, we can write for the translating observer 


P(r.t) = 1 P n e itn *- (n+q)nt)1 


(3.53) 


n =-ea 
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where 


x J n ( az <K)H^ 1) (az > K)dwdz 0 


_^ia(u-n-q)x 1 /r T 


(3.54) 


The 
rotor . 


on P n signifies that, for the moment, we are dealing with a one b laded 
The square root in the normalized radial wavenumber 


K = a|m^w 2 - (u-n-q ) 2 ( 3 . 55 ) 

has the same interpretation given above in conjunction with Figure 5 regarding 
branches of the square root. The branch points are now 


i. + _ q +n 
1 +Mx 


and 



(3.56) 


Normalization of the source function ^ brings in the effect of sweep expli- 
citly as a phase lag. We define a non-dimensional shape factor F (X) for the 
chordwise distribution of axial force per unit area on the blade as sketched at the 
top in Figure 6 . Regarding notation, note that the subscript x denotes axial 
component whereas X is the chordwise ordinate running from -1/2 at the leading 
edge to 1/2 at the trailing edge. With normalization based on local section helical 
speed, the dimensional and non-dimensional quantities are related by 


f 


X 





7 -MCA 


] 


(3.57) 


where b is blade chord. When this is substituted into the definition of the 
chordwise loading transform given by Equation 3.21 and the integration variable is 
shifted and rescaled by 


the result is 


X = ( 7 -MCA) /b 


W 


1 „ IT 2 

2 'o U o 


ik MCA 
be 7 


JFxOO 


ik bX 


dX 


(3.58) 

(3.59) 


Then, by recalling the changes of variable given by Equations 3.17 and 3.52 and 
defining the non-dimensional source transform 


**0O - fY /2 F *(X)e ikxX dX (3.60) 

we arrive at the relationship 

^[A (w ' q) ] = K U ° b e* 5 *x(k x ) (3.61) 

where the upper case ^ is used uniquely for normalized source functions. The 
phase lag due to sweep is 


_ 2 a(w-q) MCA 
a o D 

and the non-dimensional chordwise wavenumber 


(3.62) 


k 


X 


2 a(w-q) 

' a o ' B ° 


(3.63) 


where B p is the local chord to diameter ratio b/D. Recall that b is measured at 

constant radius. For 9^ and & r , the above derivation carries through with only a 
subscript change. 
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For the thickness source, the non-dimensional thickness distribution is given 
by H(X) at the bottom in Figure 6 so that the actual thickness can be written 


1.(7) - bt t h[2^) 


(3.64) 


where t„ is the section thickness to chord ratio and t b b is the maximum thick- 
ness. Substitution into Equation 3.8 and proceeding as with the loading source 
normalization leads to 


n IP 

= b Z t b e 5 * v <V 


with the normalized thickness transform 


w - til H(x)eikxX ® 


Combination of the above results with Equation 3.28 leads to 


K “ ' 0 U o «^ A 


(3.65) 


(3.66) 


(3.67) 


where 


* n ( k x> “ k x C b *v< k X ) + + B D *r< k x> £<•'> 


(3.68) 


and 




- a(w-n-q)B D ^ 3( 


— R 

D * 


(3.69) 


combines the axial and tangential force components. The partial derivative at the 
end of Equation 3.68 is in anticipation of an integration by parts when the source 
function is inserted in the radial integral. The empty parentheses (•) denote 
the, as yet undefined, material that the derivative will operate on. 

An alternative form for the combination of non-radial force components can be 
written in terms of lift and drag coefficients. Here, the reference directions for 
lift and drag are perpendicular and parallel to the local section advance direction 
as defined by the velocity triangle in Figure 4. This is analogous to the defini- 
tion used in wing lifting surface theories where lift acts normal to the flight (o 
advance) direction, and is distinct from the definition used in lifting line 
theories where the equivalent 2D lift effect is tilted back by the induced flow 
angle. Lift and drag are related to thrust^ and torque via a rotation of coordinates 
by the advance angle in the following equations 

(3.70) 


^£c 

U 


c L \ 


. V. c « 


V Or, 

** = ' U" Cl ^ l 


U 




(3.71) 


Sign conventions have been chosen so that lift and drag are positive in the usual 
sense. The s for thrust and torque are related to forces on the blades so that 
the usual thrust would result in a positive % and the usual torque would result 


Uie usual ULL1 u wuuiu — r x J 

in a negative In terms of the x wave number given above and a newly delined 

y wavenumber , the alternative form for is 


* F - \ \ C A + \ k x C d*d 


(3.72) 
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where 


K - 7^[ N-n-q)a z z o - ^]b d (3.73) 

Note that, in this representation, the integrated section lift and drag are speci- 
fied by C L and C D and the shapes of lift and drag distributions are given by f L 
and f D . These shape functions have unit area so that the ty's are defined by 

V k x> = J f L< X ) elkxX dx (3.74) 

^ k x ^ = X f D ( x ) e lkjcX dx (3.75) 

With the newly defined source function in Equations 3.67-68, the pressure harmonic 
in Equation 3.54 can be replaced by 


P 


n 



V k x ) J n ( az < K > H^ 1 } ( az>K) dwdz o 


(3.76) 


where <f> Qsx is the combined phase due to offset (Equation 3.29), sweep (Equation 
3.62), and observer axial position 


^osx = ^FA + ^3 + K (3.77 ) 

and the phase associated with the axial position of the observer is 

<j> x = a(w~u-q)x l /r T (3.78) 


Multiple Blades and Interblade Phase Angle 


In the previous section Equations 3.53 and 3.76 were derived to represent the 
pressure field of a rotor with one blade. To deal with the multiblade case, it is 
assumed that the rotor comprises B blades that are identical and equally spaced 
around the axis. In the case of steady loading, the loading is assumed to be iden- 
tical on each blade. In the unsteady case, the loading waveforms on all blades are 
identical but shifted in time. Thus, by proper choice of parameters, the loading on 
the blades of a 2-bladed rotor could be in phase or out of phase with each other. 

In general, the blade loading will be permitted arbitrary frequency, specified by q 
(not necessarily an integer multiple of the rotation frequency) and an interblade 
phase angle, specified by the circumferential order k of the distortion pattern 
(i. e., the number of nodal diameters). As will be seen, this will provide the 
generality to deal with unsteady loading caused by interference with fixed or rotat- 
Ing distortion. The fixed distortion could be caused by fuselage effects or angular 
inflow. The rotating distortion could be caused by blades of an upstream counter 
rotating rotor with a different number of blades and a different RPM. Furthermore, 
the loading could be induced by vibratory blade motion at an integer or non- integer 
multiple of the shaft rotation frequency. The extra blades are accounted for by 
simply adding their pressure fields according to linear superposition. The blades 
are identified by the counter b 2 starting at 0 for the base blade and running in 
the direction of rotation to B-l. Two effects must be included before the fields 
are added: the blades have different points of action and they have different 
phases. Equation 3.53 showed that the <^,t dependence of the pressure field for 
blade 0 is given by the exponential 

i [n<£- (n+q )nt ) ] 

e (' l iq\ 
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Moving that blade to the angular position 27rb 2 /B produces an effect equivalent to 
moving the observer to minus that angle. Thus, to account for the point of action 
effect, <f> in Equation 3.79 must be replaced by <f> - 2?rb 2 /B. For the phase effect, 
we state that the loading on blade b 2 leads that on blade 0 by an amount 27rkb 2 /B 
where 2?rk/B is the interblade phase angle and k is the number of nodal diameters in 
the unsteady loading pattern. The sign convention for the phase lead corresponds to 
the fact that, for interference with fixed distortion, a blade with the higher index 
b 2 encounters the distortion first. Thus, to account for interblade phase angle, 
we add -27rkb 2 /B to the exponential in Equation 3.53. With both effects accounted 
for, the result is 

i [n(<£- 2 ?ib 2 /B)-(n+q )ftt- 27 rkb 2 /B] ^ gQ^ 


Thus, for blade b 2 , the phase factor 

— 1 2 ?r (n~Kk ) b 9 / B ot \ 

e z (3 . ol) 

must be inserted in Equation 3.53. Adding the pressure fields of the blades results 
in 


P 2 (t) = 


B-l 




(n+q)Ht] 


-i 2 n (n+k )b 9 /B 

e z 


(3.82) 


where the subscript on p^ is to distinguish this multiblade expression from the 

earlier expression for a one bladed rotor. The subscript 2 is also needed in a 
later section where it will relate to the rear rotor of a counter-rotation system 


Equation 3.82 can be rewritten 


P 2 (t> 



i [n<£- (n+q)nt] 


B-l 

Z 

v° 


-i 2 n(n+k)b«/B 
e * 


(3.83) 


The sum on b 2 is 0 unless (n+k)/B is an integer. When (n+k)/B is an integer, say 
m, then the sum is equal to B and the double sum can be written 

p (t) = B l e 1IWB " k, *‘ ( “ (3.84) 

2 

This combination will occur so often in the acoustic and aerodynamic analysis that 
we will write 


p (t) = B 1 P n k e iCn *- (n+q)ntl (3.85) 

where it is to be understood that n=mB-k. Finally, we absorb the B into the har- 
monic definition so that the - on P n can be dropped and the complete expression 
for the pressure field in fuselage coordinates becomes 

p (t) = l P n k e i[n ^ (n+q)nt] (3.86) 

2 m=~® 

with 


^n,k 


-iaB p c 

o 


2 

o 


8tt 



i 4> n 


^ n ( k x ) J n ( aZ < K )^ 1) < aZ > K > dwdz 0 


(3.87) 


A brief discussion of the modal behavior implied by Equation 3.84 is in 
order. The exponential is 
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E = (mB - k)0 - (mB - k + q)ftt (3.88) 

which clearly represents a mode with mB-k lobes in the circumferential direction 
and a frequency (mB-k-q)0/27r as observed in the fuselage coordinate system. The 
mode is spinning about the propeller axis at a rate that can be deduced by noting 
that the zero phase point from Equation 3.88 has angular location 


<j> 

const 


mB-k+q 

mB-k 


Qt 


(3.89) 


so that the modal spin rate Is given by the same expression without the t. The mode 
can spin faster or slower than the rotor speed and in the same or opposite direc- 
tion. Specific cases for fixed and counter-rotation distortion have been discussed 
in Reference 14. 


In the blade frame, <£=<^+flt and the exponential becomes 

E = (mB-k) < 7 ^ - qftt 


(3.90) 


again indicating that qft is the only frequency observed in blade coordinates. 


Sum Over Loading Harmonics 

In this section, we bring in the effect of multiple loading harmonics and 
establish criteria for converting the double sided sums to sums with only positive 
indices. We also verify that the pressure formulas yield real results despite a 
great deal of complex arithmetic. The multiple loading harmonics are treated in the 
context of interaction with another rotor. The reason for this is that it leads to 
perfectly general results, since arbitrary loading frequencies can be generated by 
varying the speed and number of blades on the upstream rotor. For example, stopping 
the front rotor and setting its number of blades to 1 gives results applicable to 
fixed flow distortion. This report deals only with the effects of specified 
unsteady loading or specified air angle harmonics; the methodology to compute this 
input from front rotor wakes or other distortion fields has not yet been completed. 

To distinguish the 2 rotors, the front and rear will be labeled 1 and 2, 
respectively. The following equations are written for the pressure field of the 
rear rotor only because one generally thinks of wakes from the front rotor impinging 
on the rear to cause unsteady loading and noise. Of course, the front rotor also 
experiences unsteady loading in the counter-rotation case and the reader can supply 
formulas for this by switching indices. The angular speeds and blade counts of the 2 
rotors will be denoted Q 1 and Q 2 and B 1 and B 2 . When these are not subscripted, they 
apply to the rear rotor. 


In the above derivation unsteady loading occurred at only one frequency w 
For example, the loading in blade -fixed coordinates for the x force component was 

f x (7) e 1Vot (3.91) 

where u o was set equal to qQ. This can be generalized to represent the unsteady 
loading waveform as a harmonic series. The spatial frequencies in the wake from the 
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front rotor are kB 1 where k takes on all integer values. The relative speed of 
the 2 rotors is n i +fi 2 so that the frequencies sensed by the rear rotor are 
(fl 1 +n 2 )kB 1 . Thus, the waveform for the x component of loading can be written 




Z f A (7) 


k=“« 


^-i(fi 1 +0 2 )kB 1 t 


(3.92) 


We are interested in real loading only. This is assured by insisting that 
f k) - f*^ where the * denotes complex conjugate. 

From the above it can be seen that noise formulas for the general counter- 
rotation case are obtained from Equations 3.86 and Equations 3.87 by replacing q, k, 
and n as follows and summing over the loading harmonic index k. 

(3.93) 

(3.94) 

(3.95) 


q (l+Q 12 )kB 1 
k -> kB, 


n -+ mB 2 -kB 1 

where we have defined resu ^- t i- s 


P,(t) - Z Z P n.k' 

A m=-® k=-® 


i[n<£ - (n+qJHtJ 


where 


(3.96) 


-iaBp c\ 


r n,k 


8?r 


M‘ 


e W ° ,x * n t k (k x ) J | n | ( a Z< K) Hf * j (a Z> K) dwdz Q 


a). 


(3.97) 


The absolute value signs can be included or not because of the property that 
j = (-l) n J n and Y_ n = (-l) n Y n . The composite source function with k subscripts is 

» nk (V = k x t b «t v (k x ) + i* Fk (k x ) + B Ak (k x ) gf-(-) (3.98) 

' O 


with 


or 


(3 - 99) 

<3 ' 100) 


The (•) notation was explained in conjunction with Equation 3.68. Of course 
there is no sum on k for the thickness source. The wavenumbers and phase lags are 
collected here for reference: 


K = - (w-n-q) 2 (3.101) 

k x - £(„ - q)B D (3.102) 

O 

^ . i[(„.„. q)a ^ -iL] B„ (3.103) 

^.Z[ ( „. n . q)a% .^]f (3.104, 

K - tt <3 - 105) 

* - .(«•,) 5 (3.106) 

r x l T 


The notation for phase due to face alignment or offset <f> ?A was 4> q in previous 
references by Hanson. 
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Before proceeding with the discussion of folding double sided sums, we have to 
establish that 


P -n,-k = [ p n,j* (3.107) 

To this end, note that the integrand in Equation 3.97 can be written as the product 
of a real constant times 3 factors: 


[g osxj [iJ| n |(az < K)H[*j(az > K)] (3.108) 

If the signs of the n and k indices are changed everywhere and the substitution w -u 
is made, it is found that k x , k y , q, and all of the <f>* s also change sign. 

Thus, for the first factor in Equation 3.108, 


i <t>. 


osx 


~i<t> i<£ 

osx = [ e osx ] 


(3.109) 


Also, for real loading it can be shown for the source function in Equation 3.108 
that 


W -n,-k “ *n,k (3.110) 

For the third factor, recall that iJH (1) is interpreted as follows 


iJ |n |(az < |K|)H}‘f(az > |K|) for K z >0, 

-iJ| n |(az < |K|)Hj z |(az > |K|) for K 2 >0, 

| I jn|(az<|K|)Kj n |(az>|K|) for K 2 <0, 

Since H^ 2) (x) = [H^ 1J (x) ]* and J n , I n , and are 
changes sign the third factor in Equation 3.1 
From the above plus the fact that 


(n+q)>0 

(3.111) 

(n+q)<0 

(3.112) 

all (rrfq) 

(3.113) 


real, it follows that whenever (n+q) 
38 changes to its complex conjugate. 


J^f (w)dw 


/ f(-w)dw 

-CO 


(3.114) 


for any function, Equation 3.107 is established and permits us to deal with folding 
of the sums . ° 


The pressure waveform in Equation 3.96 was written in terms of a complex ampli- 
tude and an exponential. For this section only, these can be combined into one 
variable Q n k as follows 


n = P Un0-(n+q)nt] 

^n,k r n,k e 


(3.115) 


which from Equations 3.93 - 3,95 and Equation 3.115 has the property Q k =Q* 
The waveform can now be written 


p 2 (t) 

The sums can be split as follows 


Z Z 

I=-C0 v = - 


Qn.k 


(3.116) 


P 2 (t) 


< 0.0 


-uj 00 co 

+ I Qm.o + Z Qm.O + Z 


m=-l 


m=l 


Z 

k=-l 


Qn.* 


+ z 

m=-co 


Z Qn, 


k=l 


(3.117) 


where the first 3 terms deal with thickness and steady loading and the last 2 deal 
with unsteady loading. In the second and third terms we change m to -m and in the 
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third term we change k to -k with the result 

P 2 (t> - Qo.o + + <Vo> + j. + Q ».> > 

From the property just established on negative indices, this is 

P,(t) - Qo.o + 1 «£.0 + Q n ,o) + i + Qn.k) 

2 ’ m= ® k=l 


or 


P 2 <t> - Qo.0 + 2xR P - j^.o + 2XR ' P '.I. j.'V* 


( 3 . 118 ) 


( 3 . 119 ) 


( 3 . 120 ) 


where R.P. means "real part of". Thus, on returning to the original notation, the 
pressure 


CO CO 

p,(t) - l l p n 


i [n<£- (n+q)flt] 


( 3 . 121 ) 


m=— ® k=-°° 


can now be expressed in the form 


p,(t) - Po 0 + 2XR.P. £ P „ e*"”''- 01 ’ + 2XR.P. | j. P„, k e 

2 ' m=l m “ k_1 


i [n^- (n+q )Ht] 


( 3 . 122 ) 


The first term is the time -average pressure, the second is the unsteady pressure due 
to thickness and steady loading, and the third is the unsteady pressure due to 
unsteady loading. 

An alternative form can be derived by the same method 

-ikB-, t ) o n n V V u i [n</»- Cn+q)0t] 

P 2 (t) = p 0 ,o + 2 xR - p - k I p o.k e + 2xR - p -,| 1 J_ p ».* e 


/'X 10 ^ 


In either case, in the double sum, only one of the sums can be folded (eliminating 
the negative summation indices). The other must run from -« to 

In this section we have derived general expressions for the disturbance pres- 
sure in the field caused by a propeller with blade loading and thickness. This is 
in a form directly applicable to near field noise predictions, as shown in the next 
section. In later sections the formulas are specialized for far field noise and for 
aerodynamic applications. 
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APPENDIX TO SECTION 3 
SOURCE INTEGRATION IN AXIAL DIRECTION 


Part of the derivation in Section 3 required integrating a Fourier description 
of the source over the fluid volume. The x Q component of this integral was given 
by Equation 3.32, which is repeated here. 


-r 


iK, 


c M x<R+ 




n+q 

i(K x — n)x 0 


dx_ 


(3A.1) 


The integral is required for -°°<K X < 00 and -°o<rH-q<co for ^<1. This appendix shows how 
the integral is evaluated for in various ranges and deals with choosing correct 
branch cuts for the square root that appears in the radial wavenumber. 

First, we change the integration variable according to 

x o - R 0 sinh u (3A.2) 

which leads immediately to 


r® iK v M v R n cosh u + i(Kv — ri— fl)R sinh u 

I=Je xxo ^ v ° du 

-oo 

This is in the form of Heine' s formula, which is quoted from Reference 24; 

i7TH' 1) [ a|z Z -{- 2 ] = f gizcosh u + if sinh u du 

which is valid for 


(3A.3) 


(3A.4) 


lm(z±0 > 0 (3A.5) 

The requirement on the imaginary parts of z and f can be satisfied by including a 
small damping factor e in the first exponential of Equation 3A.3: 


-£ 


n+q 

i(K +ic)M R cosh u + i(K — rr-fi)R sinh u 
e -x x o x v o du 


(3A.6) 


For e > 0, Heine's formula gives 

I - i»rHg 1> (R„K r ) (3 A. 7) 

where 


K 


^(i^+ioX 2 - [^-(n+qjO ] 3 


(3A.8) 


To establish the branch cuts for K r> note that by factoring the argument of the 
square root and setting /3 = Jl-M^. It can be shown that can be written 


K r = -^OVk+XVlO 

whose zeros are 


(3A.9) 
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+ (n+q)0/V 

' " 1+Mx 


ie 


Mx 

1+M, 


(3A.10) 


and 


k' 


(n+q)Q/V 

1-Mx 


+ ie 


1-Mx 


(3A.11) 


Now, if we define 

Kjj-k* = R + e i9+ and V k ~ - R'e 10 


we can then write 

K r = /3 *Jr + R' e i(0++e 4,0/2 - |K r |e i9 


(3A.12) 


(3A.13) 


where 

„ ff + +9 +7T 

9 = 2 


(3A.14) 


The zeros k + and k~ are the branch points of the square root and the branch 
cuts can be taken as sketched in Figure 7. The consequence of choosing e > 0 is 
that k + is always below the real axis and k' is always above. Thus, if the K, 
integration is along the real axis, Equation 3A.5 is always satisfied. 

To permit integration over the source angle <f> o in Section 3, the Hankel 
function in Equation 3A.7 is expanded via the Bessel function addition theorem given 
in Reference 25 on page 363 with the result 


I 


i* 1 

m=-® 


e 


o 


(3A.15) 


where r< and r s are the smaller and larger of the source radius r Q and the field 
point radius r, respectively, in cylindrical coordinates. 

We can now safely pass to the limit of zero damping, e-0, in which case 
Equation 3.A8 becomes 


K 2 - 


[V(n+q )§] 2 


(3A.16) 


Integral I is well behaved for K* everywhere in this limit except near the 
branch points, where it has logarithmic singularities. These can be integrated 
analytically or can be circumvented by deforming the integration contour for K* 
above k + and below k‘ as shown in Figure 5. Equations 3A.15 and 16 are the results 
needed to continue the derivation in Section 3 at Equation 3.37. 
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SECTION 4 

NEAR FIELD NOISE PREDICTION FORMULAS 


In the preceding section Equations 3.96 to 3.106 and 3.122 were derived as a 
general representation of the pressure field of propeller due to thickness and 
loading. In this section the pressure equations are written in various special 
forms for noise applications. The special forms are 1) for thickness and steady 
non-radial loading, 2) for radial loading (steady and unsteady), 3) for unsteady 
loading caused by interference with wakes from any upstream counter-rotation rotor, 

4) the same where the 2 rotors have the same numbers of blades and the same RPM, and 

5) for unsteady loading cause by interference with fixed (non-rotating) inflow 
distortion. The special forms are related to equations that have appeared in previ- 
ous papers by the author. With regard to counter-rotation, the reader Is reminded 
that this report deals with noise and airloading of a single rotor caused by speci- 
fied non-uniform Inflow. The means to compute the non-uniform inflow from front 
rotor wakes or other distortion fields are not included in this report. 


Thickness and Steady Non-Radial Loading 

For steady loading, k and q are 0 so that n=mB. In this case only the second 
term from Equation 3.122 Is needed 


p (t) = 2xR . P . I e imB( *- nt) (4.1) 

m=l 

where 0 can derirved from Equation 3.87 by noting that for acoustic 

applications the field point radius ratio z will generally be larger than the 
blade tip radius so that the z <f z > notation is not necessary: 


mB 


-iaBp c? 

r o ° 

8n 


m: 


i<t>„ 


l k *W k *> + iV k *)] 


x J mB<> z o K ) ^(“K) dw dz o (4.2) 

The source transform was given by Equation 3.100: 

^F = 2 [^x C D^D^x) + ] (4.3) 

C L =C L0 and C D =C D0 are the steady lift and drag coefficients. The radial wavenumber 
from Equation 3.101 is 

K - «jM^j z -(w-mB) 2 (4.4) 

the composite phase factor, from Equation 3.77, is 

^osx = ^FA + + Ac < 4 ' 5 > 

where the phase lag due to face alignment and sweep (Equations 3.104 and 3.105) are 


<t> 


FA 



mBo Q 4 FA 
z o J D 


2 aw MCA 
S " ° 0 D 


(4.6) 

(4.7) 
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(4.8) 


The phase associated with observer position, from Equation 3.106, is 

x i 

<f> = a(u)-mB) =- 

x *- j 

The x and y wavenumbers were given in Equations 3.102 and 3.103: 


I” -a 2 z n mBcr "1 

- 2 L-*r v + sr J 


(4.9) 

(4.10) 


If the wavenumber integration variable is changed via the substitution u>=mBk, then 
the above equations become the same as those published previously in Reference 13. 


Radial Loading - Steady and Unsteady 

A formula was presented in Reference 26 for the far field noise due to steady 
radial loading without giving details of the derivation. Here we present the more 
general near field formula including both the steady and unsteady effects. The far 
field results are given in Section 5. 

The pressure waveform was given by Equation 3.96, which is repeated here for 
convenience . 


P 2 (t) -ll P n , k e i[n *- (n+q)nt] ( 4 . 11 ) 

^ m=~® k=~® 

When only the radial load source from Equation 3.98 is retained from and 
inserted In Equation 3.97, the pressure harmonic becomes 


n,k 


-iaBp 

87T 


B D * r (k x )A- [j n (az < K)H^ 1) (az > K)]dwdz o 


(4.12) 


For noise it is assumed that the field point is always beyond the tip so that z>=z 
and z<=z . Then the derivative acts only on the J Bessel function with the 
result 


r n,k 


- iaB / J 0 C o 

87T 


M 




k A( k x) ( az <K) H^ 1 5 ( a Z> K) dwdz 


(4.13) 


where 


k r - aKB D 


This result is provided for reference here. It has not been coded. 


(4.14) 


Unsteady Loading - Non-Radial 

This section presents formulas for noise caused by unsteady loading, excluding 
any radial component. Far field versions of these formulas were presented in Refer- 
ence 11 and will be given again in Section 5. 
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Three cases are treated: in the first, the 2 rotors have different speeds and 
different blade numbers. In the second, the rotors have the same speed and blade 
numbers and in the last case the front rotor is stopped so that the results apply to 
an inlet guide vane configuration or the general case of fixed flow distortion. If 
the first case alone were coded, the second and third could be computed from it 
simply by using the proper input; however, if only the special cases are needed, the 
code for them would be somewhat simpler. The following equations are written for 
the radiation from the rear rotor only but the reader can supply formulas for the 
front rotor by switching indices. The angular speed and blade count of the front 
rotor are denoted by and B^ For the rear rotor they are denoted by 
and B z or simply Q and B. 

In Section 3, it was seen that noise formulas for the general counter-rotation 
case are obtained from Equations 3.95 to 3.106 with qfl set to (Q 1 +Q 2 )kB 1 , k set to 
kB x , and n set to mB 2 -kB 1 . Results are listed below. 


General Counter-Rotation Case 


00 CO 

p,(t) = I X 

4 k=- fi o 


i[(tnB 2 -]cB 1 )^ - CmB 2 +kB 1 n 12 )n 2 t] 

P n,k e 


where the rotor speed ratio Q 12 =Q l /Q 2 and 


(4.15) 


and 


n . k 


a Vo 

87T 


m: 




* Fk ( k x> J n< az o K ) H^CazK) dv dz o 


(4.16) 

^Fk ( k x ^ = 2 [ ky^-k^Lk ( k x ) + k x^Dk^Dk ( k x) ] (4-17) 

K- ^w 2 - (w-mB 2 -kB 1 Q 12 ) 2 (4.18) 

n=mB 2 -kB 1 (4.19) 

k x - ¥l w - kB id + n i 2 )]B D (4-20) 

O 

k y = (v-mB 2 -kB i n i2 )a 2 z 0 - (mBj.-kB,)^ J B D (4.21) 

(w-mB 2 -kB 1 0 12 )a 2 z 0 - (mB.-kB,)^]^ (4.22) 

^ = |§[ W . kB.d + n i2 )]^ (4.23) 

< f > ^ = a(w - mB 2 - kB 1 H 12 )^ 1 (4.24) 


For the steady loading case, these formulas reduce exactly to those In Equa- 
tions 4.1 through 4.10 by taking k=0 only. From the form of Equation 4.15, It can 
be seen that sound harmonics appear at frequencies (mB 2 n 2 +kB i n i )/27r where m and k 
take on all integer values. 
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Special Case - Equal Blade Numbers and Equal RPM's 


In this case we can set fy-Oj-n, B 1 =B 2 =B. Then, by shifting the origin of the m 
series and substituting mBw for u, the pressure is found to be 

P (t) - Z Z p n k e i[(ro ' zkW ' mBnt] (4.25) 


and 


L n,k 


1 f® 

m 2 e^osx tf Fk (k x ) J n (mBaz 0 K)H^ 1) (mBazK) du dz Q 


z h 




(4. 

26) 

fkOS,) 

= 2[ 

( ) + (k x ) ] 

(4. 

27) 


K = Jm^ z - 

(w-l) Z 


(4. 

28) 


n=(m-2k)B 


(4. 

29) 


k - J^(mBu 

x u o 

-2kB) B d 


(4. 

30) 


l [mB(a 2 z 2 w ■ 

O O 

- a 2 ) + 2kB ] 

b d 

(4. 

31) 

TJ 

> 

11 

z 2 o [ mB(a z z 2 w 
0 0 

- a 2 ) + 2kB] 

FA 

D 

(4. 

32) 


k - ^(mBw - 
x 

2kB) 


(4. 

.33) 


<f>^ = amB (w 

X 1 

-1) r 1 

L T 


(4. 

.34) 


amB P q Cq 
87T 


The frequencies sensed by the fuselage- fixed observer in this case are given by 
mBfi, i.e. harmonics of blade passing frequency. Thus, despite the counter-rotation 
interaction, the frequencies that appear are exactly the same as for a single rota- 
tion propeller or a fan with guide vanes. The main benefit of shifting the origin 
of the m summation is that the argument of the Bessel function becomes independent 
of k. This takes advantage of the fact that standard Bessel function subroutines 
return all orders of the function for a specified argument. 


- 27 - 



Special Case - Guide Vanes or Other Fixed Inflow Distortion 


These results are obtained from the general formulas in Equations 4.15 to 4.24 
by setting B 2 =B, Q 12 =0, and B 1 =l . Also, the to integration is modified as in the 
preceding formulas . 


p 2 (t) 


00 00 

l l 

m=-« k =-co 


p gi[(mB“k)$ - mBOt] 
n , k 


(4.35) 


L n,k 


amB Z p cj; 

0 

87T 


M 




# Fk (k x ) J n (mBaz < K)H^ 1) (mBaz^) dw dz o 


^Fk^) “ 2 [ + ^x^Dk^Dk^x^] 


(4.36) 

(4.37) 


and 


K - Jm^w 2 - (w-l) z (4.38) 

n=mB-k (4.39) 

k x = |^(mBw - k) B d (4.40) 

ky = ^l[mB(a 2 z 2 W - a 2 ) + k] B,, (4.41) 

^FA - z^[>(a 2 z 2 W - "o> + k ] f (4.42) 

^ - J^mBu - k) ^ (4.43) 

4 = amB(tc>-l) ^ (4.44) 


The formulas above apply to the general fixed inflow distortion case where all 
of the loading harmonics are generated. For the guide vane case, replace k every- 
where by kV, where V is the number of inlet or exit guide vanes (not to be con- 
fused with velocity V used elsewhere in this report). Again, note that the only 
frequencies are given by mBQ, i.e. BPF and multiples. 
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SECTION 5 

FAR FIELD NOISE PREDICTION FORMULAS 


In this section the general formulas of Section 3 (Equations 3.96 to 3.101) are 
specialized, for observers far from the propeller. It Is shown that, at large ^ 
radius, the wavenumber Integral of the near field formulas can be done analytically 
by the method of stationary phase. Then the far field formulas are written out in 
various forms for the general counter rotation case and for the same special cases 
treated in Section 4 for near field noise. 

Stat ionary Phase Integral 

In Section 3, Equation 3.96 gave the pressure waveform as 

- i i p„, k (5.1) 

^ V = -eo 


The counter rotation kinematics enter explicitly via the substitutions for n, q, and 
k given in Equations 3.93 -95. The pressure harmonic is written in the following 

form 


■iaBp^c^ fl 


n * k 8?r 


J In - 


, dz„ 

k o 


(5.2) 


so that we can concentrate on the wavenumber integral I n k . This has 2 forms, 
depending on the sign on n+q: 


n ,k 


f X n,k for n+ q >0 1 
l !n,k for n+( i <0 i 


(5.3) 


These forms are (see Equations 3.111 and 3.112) 

^k - P.”- W J|n|( aZ ol K D H^}(a Z |K|) dw 

, w. 


(5.4) 


and 


i; k - - pe** 0 * e 1 " 51 * n (k x ) J| n |(az o |K|) Hj‘f(az|K|) dw 


( 2 ), 


(5.5) 


where the integration ranges specified by and w 2 encompass only the radiating 
values of the integration variable (where K 2 >0) . The offset and sweep phases 


6 = <t> + <t> 

V os r FA Y s 

were given in Equations 3.104 and 3.105 and the phase associated with observer 
position is 


(5.6) 


* a(w-n-q)x 1 /r T 

as given by Equation 3.106. The radial wavenumber (Equation 3.101) is 

K - (w-n-q) 2 


(5.7) 


(5.8) 
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We will analyze the case for n+q>0 in some detail and then simply explain the 
differences for the n+q<0 case. For large field point radius z, the Hankel func- 
tion can be replaced by its large -argument asymptotic form 


H^(y) 



e 


nn « , 

1(7 ' T ' ♦ > 


(5.9) 


With this substitution for the Hankel function, the wavenumber integral becomes 



|n|ff jt 
2 _ 4 ) 





ia[ (u-n-q)x. /r T +zK] , 

e 1 1 dw 


(5.10) 


where the branches of the complex square root of K in Equation 5.8 require special 
consideration as discussed in the appendix to Section 3. For the far field, we will 
use retarded coordinates for convenience. Thus, as suggested in Figure 8, x 
gives the position of the propeller when it emitted the signal that arrives at the 
observer when he sees the propeller at x x . The 2 variables are related by 


x i“ x r- M x'J x r + y 2 


(5.11) 


We have defined ,jx 2 +y 2 as the radiation distance. In terms of the retarded angle 
9 and the distance parameter 



R = -Jx 2 +y z /r T 


(5.12) 

the last exponential in I* k can be expressed as iRh(u) where 



h(u) = aKsin^+a(w-n-q) (cos^ -M x ) 

(5.13) 

so that 

I*. k = JgM ei ^ h<u) dw 


(5.14) 

with 

. ( . I n l w *. 




g(w) - e 1 “ 2 * ^ n (k x )J| n |(az o K)^ 

1 2 

(5.15) 


*ra(y/r T )K 


where y/r T =z . 

According to the method of stationary phase (see, for example, Jeffreys and 
Jeffreys, Methods of Mathematical Physics , Cambridge University Press), the 

integral in Equation 5.14 has the following value in the limit as R-** 


n,k 


M 


2 n 


|Rh K)l 


g( W o) 


i[Hh(w 0 )*j] 


(5.16) 


where is the stationary phase point given by h' (w o )=0 and the + sign goes as the 
sign of h" (w o ) . 

After some straightforward but tedious manipulations, the stationary phase 
point is found to be 


n+q 

W ° l-M^cos# 


(5.17) 
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in terms of the tip rotational Mach 


and the quantities needed to evaluate I* k 
number M T are 


(n+q)M T sin# 
aK ( w o> 1 -M^cos 0 

(5.18) 

h(w Q ) = (n+q)!^ 

(5.19) 

M-(l-b^cos 0) 3 

h"(«„) - - — — % - — < o 

(n+q)M^sin 2 0 

(5.20) 

so that the integral is found to be 


nr |n|n 

. A i [ (n+q) p 5 —] i<*> 

I + -1 sin# c c o 2 o ^os 

n,k a ^ ( 1 -M^cos^) 


X \( k x) J | n | [ aZ o K K)l 

(5.21) 

For the case where n+q<0, we pick up the derivation at Equation 
substitution for the Hankel function of the second kind 

5.9 with the 

. , n* »r. 

^ Z) (y) ^ # e " ~ 2 * 

(5.22) 


into Equation 5.5. The stationary phase derivation leads to the same results as in 
Equations 5.17 - 5.20 with the exception of a sign change in Equation 5.18 and a > 
sign in Equation 5.20. The result is that a combined form can be written tor both 

I's: 

Q r |n|* 

. A i [ (n+q)— - sign (n+q) r 3 i<£ 

T -1 Sin 0 * c o z P os 


!1 ' k a ^ (l-l^costf) 


- , T r In+qlz^sing "I 

x^ |n| 1-l^COS# _J 


(5.23) 


Substitution into Equation 5.2 gives the general far field pressure harmonic for all 
n+q as follows . 

nr . M*, 


-p c z Bsin0 e 

O 


iC(n+q)7 


sign(n+q ) ^- 3 


n.k 


7 i <t> 

e 05 


8^(1-1^0030) 
r | n+q | z RpSintf ”1 

W J|-i [_ npa J dz ° 


(5.24) 


where the wavenumbers are 




k * = O? ( l-t^cosfl ' q ] Bc 
_2 r a 2 Zg(n+q)M x cos^ ~1 

z o a o [_ l-l^cosfl n J 


(5.25) 


(5.26) 
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and in ^ OS = ^ FA + 0 S > the phase lags due to offset and sweep are 


FA 


2 I"" a 2 z2(n+q)M x cos5 l “I FA 

o a o l-M^os* n J IT 


(5.27) 






n+q 




l-M^cosfl 

At the end of Section 3, a discussion was devoted to establishing that 


(5.28) 


P — P 

-m,-k n,k 


(5.29) 


to verify that the equation derived for the pressure caused by a real blade loading 
is real and so that formulas for folding one of the sums in Equation 5.1 could be 
established. Here, we want to verify that this property has been correctly pre- 
served in the passage to the far field. To accomplish this, first note from Equa- 
tions 3.93 - 95 that, when m and k change sign, so do n and q. Also, in 
Section 3 it was shown that Then, since <f> Qs changes sign with n and q, it 

is clear that Equation 5.29 holds. Furthermore, it can be seen when P n k is substi- 
tuted into Equation 5.1 for the radiated pressure field, that the exponential 

i(n+q)^p-(r - c Q t) 

e ° (5.30) 


appears for all n and q. This is the required behavior representing outgoing 
waves for every harmonic mode. 

Equation 5.24 is the general expression for the far field pressure harmonic in 
terms of the general mode index n and the frequency ratio q. To predict wave- 
forms, P n k can be used either in Equation 5.1 or in Equation 3.122, which is 
repeated here for convenience. 


P (t)- P 0 0 + 2xR.P. Jp n 0 + 2xR. P . ? 


m=l 


l 

k-1 


i [n<£- (n+q)0t ] 


n.k 


(5.31) 


In the next sections, radiation formulas are displayed explicitly in terms of rotor 
blade numbers and speeds via the substitutions in Equations 3.93 - 95. 
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Thickness and Steady Non-Radial Loading 

For this case k and q are 0 and we need only the second term from Equation 

5.31: 


P 2 (t)= 


2xR . P . 


03 





where the harmonic P m = P n 0 is 


(5.32) 


-p c^Bsin# 



x 


S/r^l-M^os# ) 

-i 

I ^(k,) 



mBz Q M T sin^ 

l-M^osfl 



(5.33) 


and the source term is 


^mB + 2 [ + ] 


The x and y wavenumbers are 


k_ = ^ , .. mB , B n 


y z, 


° 0 l-M^os# D 


2mB T ”1 p 

z 0 <7 0 |_ l-f^cosfl J D 


and the phase lags due to sweep 


and offset in <f> =<£-+<£ are 

r os FA r s 


/ = 2a mB MCA 

cr Q l-M x cos^ D 

2mB r 1-g^cosg “I fa 

FA Z 0 (7 q I 1-f^COS^ J D 


(5.34) 


(5.35) 

(5.36) 


(5.37) 

(5.38) 


These are equivalent to the far field formulas originally presented in References 11 
and 12. Note, however, that there has been a change in sign in the definition of 
k y . This is a change in sign convention only, not an error. In the references 

just mentioned there was considerable discussion of the effect of geometry, includ- 
ing particularly sweep, on noise. 
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Radial Loading - Steady and Unsteady 


Theory for this loading term was presented by the author in Reference 26 for 
steady loading only. The general unsteady loading harmonic can be written down 
immediately from Equation 5.1 and 3.98 in which the radial source component of 





(5.39) 


which acts on the Bessel function, symbolized by the empty parentheses. When this 
is inserted into Equation 5.24, the pressure harmonic is found to be 

Hr Inin 

i [ (n+q )- sign(n+q ] 

-p c~Bsin#e ° 

p = 

r n , k V 

Stt^CI-M^cos#) 


n 2 r e 03 


VrkCKc) 


r l n +q losing ~j Tt r 
[_ l-M^ostf l n l 


| n+q | z o M r sin0 
l-M^cosfl 


] 


dz 


(5.40) 


This has never been coded because unsteadiness in the radial loading source has not 
been evaluated. However, for the steady effect, the result can be written down 
immediately: 


P_ - 


imB [~ -] 

-p c^Bsin# e ° 


8^(1-^0080) 

| l%e*” M r k x sin^ (1/2) ^ r (k x )J^ 




I^cos# 


] 


dz„ 


(5.41) 


Up to this point in the derivation, has represented a radial force distrib- 
uted over the blade surface. If now we say that the radial loading the source acts 
only at the tip, then can be replaced as follows 

*r - (5.42) 


C R is a coefficient of radial force defined on the next page. When Equation 5.42 
is substituted into Equation 5.41, the result is 

imB[- — -) 

-p c^Bsintfe ° 

p - — 2_? 

m y 

8tt^( 1 -f^cos#) 


7 1 

ne 


M r k x sin0 


r mBMpSin^ 


f *.<**> J - 


] 


(5.43) 


where k x and ky and the <{>' s are the same as in the preceding section. This 

formula was presented in Reference 26 in conjunction with discussions of the radial 
loading associated with formation of a tip vortex. 

In Equation 5.43, the interpretations of C R and >l r are as follows. Say 
that s(7 o ) is the radial force per unit length distributed along the tip chord 
line. Chordwise distance is measured approximately here along the helical advance 
coordinate which was used in Section 3. The radial force coefficient C R is 
defined according to 
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s(7„) “ (P a Ui/2)r T S[( 7 Q -MCA)/b] C R 


(5.44) 


where MCA is the mid chord alignment, or sweep, at the tip. The function S gives 
the shape of the chordwise distribution and is normalized to unit area according to 

/ , S (7 /b) d( 7 o /b) =1 (5.45) 

chord 0 0 


Thus, the integrated radial force is (p U*/2)br T C R where U T is the tip helical 
speed. 
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Unsteady Loading - Non-Radial 


In this section, formulas are presented for calculation of noise caused by 
interaction of a rotor with the unsteady flow field of another rotor or with a fixed 
(non-rotating) distortion pattern. The 3 cases treated are the same as those in 
Section 4 for near- field noise. These are given by Equation 5.1 and 5.24 with q, k, 
and n specified by Equations 3.93 - 95 and the non-radial source term k i^ Fk . As 
in Section 4, the reader is reminded that this report deals with the counter- 
rotation problem only up to the point of computing unsteady loading and noise caused 
by specified non-uniform inflow. The means to compute the non-uniform inflow from 
wakes of an upstream rotor or other flow field distortion are not included herein. 

General Counter-Rotation Case 


p (t) = Yj X p - CmB 2 +kB i n i2 )n 2 t] 


k=-® 

where the rotor speed ratio 0 12 =0 1 /r2 2 and 


fl IT 

i I - signCmB^kB^fl^) |mB 2 "kB 1 |n/2] 


- ip c?B«sin0e 

r n O 2 


n .k 


Stt^I-H^cos^ ) 
x * Fk (k x ) J |mVkBi 


I” l mB 2 +kB l Q 12l Z o M T slng "I . 

L 1-M^cos# J Z ° 


- l[ k y C L k \ t (K) + k x C oAk (k x )] 


(5.46) 


(5.47) 

(5.48) 


and 



mB 2 -fkBiQj 2 

l-^cosf? 


kBi(i+n 12 ) b, 


] 


i -2 V a2z o( II ^2 + ^ B l fi 12) W x COS ^ 

^ " V'o |_ 1-MjjCOS^ 

k 2 r 

^FA = 2^1 [_ l-M^COStf 


(mB 2 -kB 1 ) J B d 
( mB 2 -kB 1 ) 


(5.49) 

(5.50) 


(5.51) 


K 



mB 2 +kB i n i2 

1-H^cos^ 


kB, (1+Q 12 ) 


] 


MCA 

D 


(5.52) 


From the form of Equation 5.46, it can be seen that sound harmonics appear at 
frequencies (mB 2 0 2 +kB i n i )/27r where m and k take on all integer values. For the steady 
loading case, these formulas reduce exactly to those in Equations 5.32 through 5.38 
by taking k=0 only. Also, by the use of the formulas J_ n (z)=( - l) n J n (z ) and 
J n (-z) = (-l) n J n (z) , Equation 5.47 can be shown to be equivalent to the forms given in 
various AIAA papers by Hanson (eg. reference 14), which do not employ the absolute 
value signs on the mode order and Bessel function arguments. 
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Special Case - Equal Blade Numbers and Equal RPM's 


In this case we can set n^-n, B 1 =B 2 =B . Then, by shifting the origin of the m 
series and substituting mBu for w, the pressure is found to be 


P 2 (t ) 





k=-co 


p i[(m-2kW - mBflt] 
r n,k e 


(5.53) 


i --sign(m) |m-2k |nB/2] 

-Ip c 2 Bsin0e 0 

r o 0 

87^(1-^0030) 

ImlBz^sin^ “j ^ 
l-t^cos# _| ^ Z ° 

(5.54) 


x| M^e os ^ Fk (h x ) d |m- 2 k |b ^ 

z h 


* Fk (V “ + k x c DkV(k x ) ] 


(5.55) 


and 


, 2a f mB 
X " 


l-M^OSI? 


2kB 
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.-ir “ B( °^ C0S ^H + 2kB 1 B d 

Z o°o L l-t^COS# J D 
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~i ia 

j ° 


(5.56) 


(5.57) 


(5.58) 


i _ 2a f niB _ ou 


^ MCA 

J D 


(5.59) 


The frequencies sensed by the fuselage -fixed observer in this case are given by 
mBO, i.e. harmonics of blade passing frequency. Thus, despite the counter-rotation 
interaction, the frequencies that appear are exactly the same as for a single rota- 
tion propeller or a fan with guide vanes. The main benefit of shifting the origin 
of the m summation is that the argument of the Bessel function becomes independent 
of k. This takes advantage of the fact that standard Bessel function subroutines 
return all orders of the function for a specified argument. 
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Special Case - Guide Vanes and Other Fixed Inflow Distortion 


These results are obtained from the general formulas in Equations 5.46 to 5.52 
by setting B 2 =B, fi 12 =0, and Bj-1. Thus, n=mB-k, q=k, and n+q=mB. 
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The formulas above apply to the general fixed inflow distortion case where all 
of the loading harmonics are generated. For the guide vane case, replace k every- 
where by kV, where V is the number of inlet or exit guide vanes (not to be con- 
fused with V used for velocity elsewhere in this report). Again, note that the 
only frequencies are given by mBO, i.e. BPF and multiples. 
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SECTION 6 

SOUND POWER AND WAVE DRAG 


Sound power radiated from a propeller might be of interest as an overall mea- 
sure of the effect of geometry changes on noise. However, because of the subjective 
nature of noise, this is a somewhat academic application. A more practical applica- 
tion is in the area of performance. Sound power is related to wave drag. Consider 
the sketch in Figure 9 showing a propeller in flight surrounded by a cylindrical 
surface S. The surface S comprising 2 discs, S t and S 3 , and the curved surface, S 2 
the curved surface, S 2 is reminiscent of figures used in aerodynamics text books 
for pressure drag analysis (e.g. Ashley and Landahl (ref. 27, page 174). The momen 
turn flux through S 2 corresponds to radiating energy and wave drag; the flux 
through S 3 is related to energy in the trailing vortex system and is used to 
compute vortex drag. (Other names for vortex drag are drag-due -to -lift and induced 
drag.) In this section, we derive the formula for sound power (or wave power) 
passing through S 2 and compare it to the shaft power. We also examine the spec- 
trum of sound power and find that it is smeared in frequency bands centered on the 
blade passing harmonics, as would be expected from Doppler considerations. The 
width of the bands increases with flight Mach number. 

The derivation herein is an extension of the analysis in chapter 11 of Morse 
and Ingard's book (ref. 28) which deals with power from oscillating point dipoles in 
straight line motion. In the present analysis, the cylindrical surface S 2 of 
radius r is constructed with its axis coincident with the propeller axis and flight 
path. The sound power II is the integral over the surface of the radial component 
of the acoustic energy flux averaged over time: 

n = iT -T* t-T pUr dt rd ^ dx 

where u is the radial velocity component which will be obtained from the acous- 
tic pressure p via the momentum equation. The averaging time T need not be speci- 
fied, since it will be found to drop out of the analysis. An important feature of 
the Morse and Ingard analysis is that Equation 6.1 is written in the fluid- fixed 
coordinate system, removing ambiguities in the interpretation of sound power that 
could be introduced by using moving coordinates. This approach is retained in the 
present analysis. 

To evaluate Equation 6.1, we obtain a suitable expression for p , then derive 
an expression for u r from it, and finally perform the <j > , x, and t integra- 
tions analytically. The expression for pressure, obtained from Equations 3.96 with 
steady terms only and 4.2 is, for the observer translating with the propeller, 


P(t) 
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The wavenumbers are 


2 aw Q 

~°7 B ° 


K = z^(' a2z l w + b d 

and in <j> os = <l> Ff +<l> s the phase lags due to offset and sweep are 
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The Hankel function 
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( 6 . 6 ) 

(6.7) 

(6.8) 

(6.9) 

( 6 . 10 ) 


To return to the fluid-fixed coordinate system, we simply substitute x x -x-Vt as 
was discussed in conjunction with Equation 3.46. In the new expression for pres- 
sure, we also split the upper and lower halves of the series so that the Hankel 
functions can be displayed explicitly. 
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where the wavenumber integration is restricted to the range corresponding to radiat- 
ing waves via 
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The second series can be rewritten by substituting m for -m and w for -w with the 
result 
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(6.13) 


The negative Bessel function orders are eliminated because Bessel functions 
and Hankel functions are odd in order. Also, the source transform with negative 
argument and index results in the complex conjugate, as shown. 


To find the corresponding form for the radial component of velocity, we inte 
grate the linearized momentum equation 
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with the result 
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Since z appears only in the Hankel function and t appears only in an exponential, 
the operations are easily performed on Equation 6.13 to obtain the working form for 
the radial velocity: 
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(6.16) 

We are now prepared to perform the integration in Equation 6.1 using Equation 
6.13 for p and Equation 6.16 for u r . To keep the m's and w's separate, we 
place primes on m and w in the equation for u r . We also adopt a special 
notation of P's and U's that combine the constants plus all of the terms under each 
integral except for the exponentials. With this shorthand, the power can be written 
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The integrand is the product of 2 series in the first square brackets and 2 more 
series in the second square brackets. When these are multiplied, 4 terms result. 

For the <f> integration, the product of the first series from each pair produces 

I 2 * e i(m+n’)* d J (6.18) 

J 0 

which is 0 since both m and m' take on only positive values. The same remark 
applies to the product of the second series from each pair. However, for the cross 
term , we have 
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This is 2 tt for m'=m and 0 otherwise. With this integration completed and the 
summations simplified, we find 
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This can be rearranged by moving the w integrals outside the z Q integrals and 
moving the x integration inside: 
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The x integrations each yield (2nr T /a) 8 (w-v' ) so that the w' integration becomes 
trivial. Furthermore, the integrand is then time independent so that the time 
integration (1/T) /dt = 1 and we now have 
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We can now return to the original notation, eliminating the P's and U's: 
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At the analogous point in the Morse and Ingard derivation for translating point 
dipoles, asymptotic forms for the Hankel functions were used to pass to the far 
field. However, this is not necessary if we note that the Hankel functions appear 
in a combination that can be evaluated using the Wronskian relationship 

wtH^Cx), H^ 2 ) (x) ] - H^>(x) A l< 2) (x) - H^ 2 ) (x) A h^>(x) 

" # (6-24) 
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The final result for sound power is 


n - irjl;' Mje‘*”»,*<v=,>J.,(az 0 K><Jz 0 | ! dv, (6.25) 

4,rM x m=l “l 1 Z h 

After a somewhat complex derivation, this is a remarkably simple result. The for- 
mula is easily programmed and can be evaluated quickly even on personal computers, 

if desired. 

Various reference powers or normalization schemes could be used for the sound 
power in Equation 6.25. However, since the prime motivation for investigating sound 
power has been to understand its relation to performance, it makes sense to use the 
same reference as in the definition of propeller (shaft) power coefficient: 

r _ SHP (6.26) 


Thus, the coefficient of sound power becomes 

r n (6.27) 

b SP “ ^ xt 3t,5 


With this definition, we have 
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This can be compared with shaft power coefficients, which are typically of order 1 

Equation 6.25 is in an undesirable form for calculations for low flight Mach 
number, i.e. as M x -0. A simple change of variables can be used to provide a 
suitable form and, at the same time, serve as a partial verification of the power 
formula. If the substitution 


= mB 
= l-M x cos^ 


is made, we find that 
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which approaches mBM^infl as M^O. Thus, the Bessel function argument 

az Q K -► mBz^sinfl 
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For small M^ the differential 


dw -► -mBM^infl dd (6.32) 

in which the M x cancels the 1/M^ in the constant of Equation 6.25. The integra- 
tion limits w 1 and w 2 become n and 0. Thus, with the substitution for w , the 
expression for sound power at low Mach number becomes 
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(6.33) 


Mach number cancels from the expression at low Mach number, and the integrals are 
obviously finite. 


To check the above result, we have used the expression from elementary acous- 
tics books for sound power of stationary sources 


n = flViT/I P 2 rZ sin8 de ^ (6.34) 
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where p is the far- field pressure and the overbar denotes the time average. The 
reader can verify that, when Equation 5.33 with M x =0 is substitued for p, Equa- 
tion 6.33 is the result. This is not a complete verification of the general results 
of this section, but it does provide some confidence that the algebra was done 
correctly. 
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SECTION 7 

DOWNWASH EQUATIONS FOR AIRLOAD CALCULATION 


This is the first of several sections on theory for aerodymanic blade loading 
and propeller performance. The main task of this section is to develop an equation 
giving steady or unsteady downwash at the blade lifting surfaces as an integral of 
Slbfade loading in the form of a standard lifting surface kernel function equa- 

tion: 


| - a(x,y) - JjAC p (x o ,y o )K(x,y;x 0 ,y o )dx 0 dy o 


(7.1) 


V 

Here AC is the coefficient of lift pressure, which is to be found, and a is the 
downwash angle, which must match the known angle between the blade camber surface 
and the blade section advance direction. In the lifting surface context, lift is 
that force component acting normal to the local blade section advance direction i 
analogy with wing methods. This is distinct from the system used in } l 

methods where thf equivalent 2-dimensional lift vector is tilted back by the induced 
flow angle. K is the kernel function giving the downwash per unit loading It has 
various forms depending on the problem being solved. For subsonic wings the we 

known kernel is 
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(7.2) 


x and y are the streamwise and transverse source coordinates in the wing surface and 
x° and y are the corresponding field point coordinates. Equation 7.2 and its super- 
sonic counterpart are presented in Reference 4 with a discussion of application 
methods. The propeller kernel is more complicated but is analogous in that 2 terms 
appear, one of which is independent of the streamwise coordinate and the other of 
which contains the streamwise coordinate and the compressibility effect ( ere via 
the 8 ) . Later sections show how these integral equations can be inverted to de 
mine airloading as a function of the blade boundary conditions. 


Acceleration Potential Method 


The derivation of this 
given in standard textbooks 
is basically an integration 


section is an extension of the pressure potential method 
such as Bisplinghoff , Ashley, and Halfman (ref. 3). It 
of the momentum equation 


D v 
Dt 


- J(F - Vp) 


(7.3) 


which relates the disturbance pressure, body force, and velocity^ DY/ Dt 
acceleration following a fluid particle. F is the body force, which doesn t ordi- 
narily enter in wing theory because the air particles pass over or under the wing. 
For propellers, the particles pass through the disc so that the role of t e 0 y 
force might seem different in this case. However, from the sketch in Figure 10 i 
can be seen that a particle whose acceleration is being integrated also passes over 
or under the blade, not through it. (The particle motion shown m the figure is 
shown as undeflected by the blade loading for simplicity. This undeflected path 
also used for the integration described below.) Therefore, the body force is not 
needed tf on^y the downwash is to be calculated. As with wing theory, the downwash 
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equations are not valid on the wake sheet, but can be evaluated in the limit as the 
field point approaches the sheet. For a complete treatment of the propeller theory, 
the body force is useful in finding the axial induction and relating it to classical 
propeller momentum theory. It also helps in interpreting the above-mentioned limit- 
ing process and therefore is included in the analysis. 

Our use of Equation 7.3 involves linearization at various points in the deriva- 
tion. First, we will use an expression for the disturbance pressure derived in 
Section 3 based on the linear wave equation. Second, the density p is set equal to 
its ambient value p°. This is justified because F and p are first order distur- 
bance quantities. Perturbations of p would Introduce second order quantities in 
the momentum equation. This is discussed more formally in Reference 3 . Third, in 
the integration of Equation 7.3, the particle is followed along undisturbed stream- 
lines as in wing methods and explained in the next paragraph. This is a good 
approximation for typical cruise conditions where the axial induction is small 
compared to the flight speed. However, for takeoff and climb conditions, the 
approximation is not acceptable and a non-linear modification must be made as 
explained in Section 8 . 

w I nte g ra t;ion of Equation 7.3 could proceed by using the expression for distur- 
bance pressure given by Equations 3.96 and 3.97, taking its gradient, and performing 
the integral. However, it is more convenient to integrate the pressure and take its 
gradient afterwards. We write the functional dependence of the pressure as 
P(x,r,^,t) in a coordinate system translating but not rotating with the propeller. 

(In this and following sections where we are dealing with propeller aerodynamics, as 
opposed to acoustics, the sign convention for x will be changed from that in the’ 
previous sections to positive downstream.) The quantity 

lKx,r,*,t) = jr P(x,r,*,t) (7.4) 

r O 

is commonly called the acceleration or pressure potential. If we write $ as the 
disturbance velocity potential, then 

$(x,r,^,t) = ^(x,r,^,r) ( 7 . 5 ) 

because taking the gradient of both sides returns Equation 7.3. The integration is 
from upstream infinity (-®,r,^) to the field point (x,r,^) occupied by the par- 
ticle at time t. A dummy time variable r is used for the integration to follow the 
particle. Thus, at any r prior to t, the particle was located at [x-V(t-r) ,r, 0 ] . 

The integral of Equation 7,5 is ~ 


$(x,r,^,t) = J^(x-V(t-r),r,*,r) dr ( 7 . 6) 

The integration variable is changed from r to x'=x-V(t-r) so that 

*(x,r,*,t) = ifjKx'.r.^.t - ^-)dx' (7.7) 

Bringing the V out from under the integral is a significant linearization that will 
be modified in Section 8 . Evaluation of Equation 7.7 and the associated downwash 
velocity comprises the remainder of this section. 
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Velocity Potential 


From now on we will deal with only one loading harmonic at a time. The 
unsteady loading frequency is qQ/27r and the harmonic loading index k gives the 
interblade phase angle 27rk/B as explained in conjunction with Equation 3.80. For 
the harmonic with index k, Equations 3.96 and 3.97 can be combined to give 
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We are analyzing only the effect of the non- radial loading source so that the source 
transform 'J^=i'I r p . Recall that the radial wavenumber is 

K = ^w 2 - (u-n-q ) 2 (7.9) 

where the branches of the square root are explained in conjunction with Equation 
3.36. 4> is the phase lag associated with the offset or Face Alignment of the blade 
section at radius ratio z Q and <f > ^ is the phase associated with the axial position of 
the field point, which was given by Equation 3.106. However, with the new sign 
convention for axial position, this becomes 

^ = -a(u-n-q)^- (7.10) 

x 1 x 

In the harmonic summation, the effects of multiple blades and arbitrary interblade 
phase angle are accounted for by summing on m but recognizing that the mode order is 

n=mB-k (7.11) 


as explained in Section 3. 

At this point it is desirable to eliminate the source transform and return to a 
direct representation of the lift pressure in terms of the coefficient AC p ( 7 o ,z 0 ), 
where 7 is the coordinate measured along the section advance direction as shown in 
Figure 4. Within the linearization approximation, 7 q is also the chordwise distance 
variable. Since we are dealing with only the lift component of the blade force at 
this point, Equations 3.72 and 3.74 show that 

MM - jVJA (x)e ‘ Vdx (712) 


where the chordwise wavenumber is 
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(7.13) 


X is the non-dimensional chordwise distance measured from mid- chord which can be 
returned to the 7 q system by Equation 3.58 to give 
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If we define a new non-dimensional distance based on tip radius 


V r T 


(7.15) 


we find, recalling C Lk f L =AC , that 
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We have defined a new wavenumber corresponding to k^ 
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and for future use we 

also define 
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Substitution of Equation 7.16 into 7.8 gives the final form for the pressure har- 
monic 
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(7.21) 


Now in applying Equations 7.4 and 7.7 to find the velocity potential, note that 
the operations are only on the x and t variables in Equation 7.21. Since these are 
only in exponentials, we can extract them temporarily: 
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With the change of variables required by Equation 7.7, this becomes 

i [ - (n+q)tlt+(n+q)^-au|— I 

E= e V 1 


(7.23) 


leaving us with the improper integral 


I - J* e rT dx' (7.24) 

-cO 

The standard method for evaluating this kind of integral is to add some damping in 
the exponent to suppress the oscillations at evaluate the integral analy- 
tically, and then take the limit as the damping goes to zero. The result is given 
by Champeney (ref. 29), for example. In the present notation, 
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(7.25) 


_ i aux 

1 = ^[ 7rfi ( w ) + ^] e T 

Thus, the effect of the integration in Equation 7.7 is simply to insert the factor 
i[jr5(w)+|] under the w integral in Equation 7.21. The velocity potential is 
obtained by including the factor -1/p from Equation 7.4 and performing the u 
integration that is facilitated by S(w). The result is the sum of 2 terms as 
follows . 
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(7.26) 


where we have introduced the shorthand notation 

JH^u) - J n (az < K)H^ 1 ) (az > K) (7.27) 

The frequency factor e" lqnt , that is common to both the loading and downwash, has been 
dropped from the expression for velocity potential in the remaining text, and is now 
implicit. 


Downwash Velocity 


The desired downwash velocity w/U is obtained by differentiating $ in a direc- 
tion normal to the helicoidal sheet according to 
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The partial derivatives for the chain rule can be evaluated from Equations 3.11a and 
3.12a with the result 
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The operations are straightforward and produce 
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^ is the angle between the advance helicoid of the load point and the advance 
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helicoid of the field point: 


2 A . Ox nt . 2a o FA . , Ox 2a o FA 

and Ax is the axial distance from the load point to the field point x. 
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Equation 7.30 can now be written as an integral equation in the form of Equation 
7.1: 
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where the kernel function is 
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Equations 7.33 and 7.44 are the downwash equations in their most general form 
giving the induced flow angle a as an integral of the lift pressure coefficient 
AC p . For steady flow, q=k=0, the 2 terms of Equation 7.34 correspond exactly to 

the 2 terms of the wing kernel in Equation 7.2. The first term is independent of 
streamwise position x and is incompressible, i.e. the speed of sound doesn't appear 
via K* or ft. The second term gives the x dependence and has the compressibility 
information via M x in the definition of K in Equations 7.9 and 7,27. Equation 

7.34 is valid at any field point $,r,x except on the advance surfaces behind the 
point where the loading begins (normally the leading edge) . On the advance sur- 
faces , the kernel has the same singularities as the wing kernel of Equation 7.2. 
These singularities are not apparent from the form of Equation 7.34 but arise from 
the non- convergence of the m-series. Dealing with the singularities is the main 
challenge in applying the propeller kernel. For accurate results, they must be 
exposed and treated rigorously. Methods for dealing with the leading singularity 
are described next . 


Wake and Bound Components of the Kernel Function 

The 2 terms that have appeared in both the wing kernel and the propeller kernel 
include a second order singularity in the spanwise coordinate. That is, when the 
field point (y in Equation 7.2, z in Equation 7.34) approaches the load point (y 

J O 

2 

or z o ) , the kernel diverges as l/(y-y ) . This behavior is actually only a property 

of the wake, as will be shown, and thus for our purposes, a better subdivision of 
the kernel is into a wake term and a bound term rather than the current form. For 
clarity, this approach will be described first in terms of the wing kernel and then 
applied to the propeller kernel. 

In examining the propeller kernel, the load is placed at the origin so that 
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Equation 7.2 simplifies to 


K 


8*y 2 £ 1 


i 2 + p 2 y 2 


] 


(7.35) 


The behavior of this is, of course, quite different far upstream and far downstream. 
For x»0, the second term in the brackets approaches 1 and K-* l/47ry 2 . For x«0, 
the 2 terms cancel and K->0. For finite x, the term in brackets in Equation 7.35 
can be expanded as a binomial series with the result that 


K -+ 


1 

4?ry 2 


1 

16 tt 



+ terms of order p 2 y z /x 2 


x>0 


(7.36) 


K -» 



x<0 


Thus, the singular behavior occurs only downstream of the load, as expected from 
physical considerations. The 1/y 2 singularity can be confined to a wake term by 
adding sign(x) to the first term in brackets and subtracting it from the second. 
Since 1 + sign(x) is twice the unit step function H(x) , the kernel now divides into 

a wake term K w and a bound term K b 


K - K w + K b 


where 


K = 


1 

8xy 2 


K 


[ 


4xy 2 


H(x) 


J7 


2 + p 2 y 2 


- sign(x) 


] 


(7.37) 


(7.38) 

(7.39) 


For large |x| , the bound term approaches -/S z /(16xx 2 )sign(x) and thus trails no wake, 
as its name implies. Furthermore, the wake term contains none of the second order 
singularity. 


To use this same method to divide the propeller kernel into a wake term and a 
bound term, note in Equation 7.34 that the integrand in the second term contains the 
factor 

F(w) - i g-i.au Ax (7.40) 

where Ax, the non-dimensional axial distance between source point and field point, 
was given by Equation 7.32. Function F controls the behavior of the w integral for 
large x. Also note that the integral of F over w is proportional to 

r“ 

sin^awAx dw _ J sign (Ax ) (7.41) 

J o 
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giving us a means to eliminate the wake behavior from the second term in K. 
Specifically, we subtract K f0 K f JH^O) F(w) from the integrand in the second term of 
Equation 7.34 and add its integral to the first term. The integral produces 
sign(Ax) with a coefficient exactly equal to the first term in Equation 7.34. This 
yields a unit step function in the same manner as the wing kernel so that the kernel 
for the propeller can be written 


— — w — b 

K = K + K 


(7.42) 


where 


(7.43) 


(7.44) 

This is the final form of the propeller kernel before discretizing in Section 
9. Like the wing kernel in Equations 7.38 and 7.39, the wake term is incompres- 
sible. Also, it can be seen that far downstream, K* is independent of x except for 
a convected wave in the unsteady case. The bound term contains the compressibility 
effects, as mentioned above. 


Bcr 


K ’ 4^T H(AX) J/ n * n eiaq ' X V? ^ n ( l n+ q|az < )K( |n+q|a Z> ) 


and 


— b 

K 


-B <7 


167TCr 2 ZZ„ 


Y e in4> e iaqAx l 


iauAx 


[K^K { JHJw) - K^K f JH„(0) ] dw 


Approach of Field Point to the Vortex Sheets 

We still need to deal with the singular behavior of the kernel as the field 
point approaches the advance helicoid downstream of the load points. To simplify 
the discussion, we do this in the context of steady loading, i.e. for q=k=0, how- 
ever, the same process applies to the unsteady load case. Also, it is instructive 
to determine the result that would have been established if the body force had been 
included in the integration of the momentum equation. This helps interpret the 
singular behavior of the kernel and also is needed in relating the circumferentially 
averaged results to classical propeller momentum theory. 

To rewrite Equations 7.43 and 7.44, we eliminate the K's defined in Equations 
7.19 and 7.20 so as to expose the summation index explicitly. The steady loading 
kernel terms then become 


-SW 
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Bcr;* 
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47TZZ^ 
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H(Ax) 


l n 2 e inJ I n (|n|a Z< )K(|n|a Z> ) 

m=-co 


(7.45) 


— sb 
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167TCT 2 ZZ 


l e 
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in# 


[K ?0 K ? JH^w) 


- K f0 K ? 


JH„(0) ] 


- i at>Ax 


dtj 


(7.46) 


where the s in the K superscript denotes steady loading. These expressions are 
further modified by taking advantage of symmetry in n in the wake kernel and by 
separating the n=0 term in the bound kernel . 


- 52 - 


(7.47) 


— sw Bcj 3 


K = 


2-ttzz 


H(Ax) £n 2 cos(4) I n (naz < )K(naz > ) 


m=l 


sb — sb — sb 

-7 

v 0 


K = K n + K n0 


(7.48) 


where 


and 


K f = Ba °* wI 0 (aj8uz < )K 0 (a^z > )sin(awAx)dw 

» n 


(7.49) 
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n0 167TC7 2 ZZ 0 


[K 40 K. JH,» - K ?0 K e JH„(0) ] s — 


p, - i awAx 

e dw 


(7.50) 


The prime on signifies that the m=0 term is not included in the sum. 

In the last 4 equations, the summations are obviously harmonic series in the 

— sb 

circumferential angle $ with harmonic order n=mB. In the bound kernel, K 0 , as 
the zeroth term in the series, represents the circumferential average induction and 

~SW , 

K Sb represents blade -to -blade variation. In the wake term K , however, we note that 
there is no zeroth element. This appears to suggest that there is no mean swirl or 
axial induction in the wake, which is obviously incorrect. Formal evaluation of the 
circumferential average of the wake term would require the integration: 


C - iC r " # <7 ' 51) 

The key to the apparent paradox is that the wake kernel is singular on the wake 
surfaces, and without the correct limiting process, cannot be integrated through the 
full range 0 to 2rr in the obvious way . 

Before proceeding, we will verify that the behavior of the wake term for 
S=0 in Equation 7.47 is correct in comparison with the wing kernel behavior in 
Equation 7.38. This is accomplished by replacing I n and K„ in Equation 7.47 by 
the first terms in their large -argument asymptotic forms from Reference 25. 
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*j27rn (l+a 2 z 2 ) 
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(7.52) 


A = *Jl+a 2 z 2 + In az - 

l+*]l+a 2 z 2 


This leads to the sum 
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which is a standard form yielding 
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(7.53) 


(7.54) 


(7.55) 


- 53 - 


When this is expanded about z-z Q , the leading terra is found to be H(Ax)/[4tt(z-z 0 ) z ] , 
This is exactly the same as the leading term in the wing kernel singularity (Equa- 
tion 7.36), as might be expected from physical considerations, and can be considered 
a limited verification of Equation 7.47. 

For the integral of the second order singularity, just exposed, we must take 

the Mangier principal value as discussed by Ashley and Landahl 27 . The method we 
will use is different from most of these and was partially inspired by an analysis 
of Goodman (Reference 30) of marine propeller wakes. We return to the integral 
equation (7.33) and insert the wake kernel to get the wake term of the downwash 
associated with steady loading 

a sw = //AC p ( 7 o ,z o ) K SW dr o dz o (7.56) 


B 

2nz 


UJ 

s ( z 0 ) I nZ c °s(n$) I n (naz < )K n (naz > )dz 0 


m=l 


(7.57) 


where 


S(z 0 > " RAC p ( 7o ,z o )H(Ax)dr o (7.58) 

is the radial distribution of the loading source. We now eliminate the z <} z > 
notation and split the z Q integration into its upper and lower portions. If we 
consider only the contributions of loading in a band of radii from z x to z 2 
spanning the field point z, the downwash is 


Q 12 W= 2%z ! nZcos ( n ^) 


f 2 

K^naz)] S(z D )I n (naz o ) dz o 

Z 1 

+ I n (naz)| S(z o )K n (naz o ) dz o I 


(7.59) 


Next, we recall that the differential equation defining modified Bessel functions 
can be written 


w(z) 


d/ z dW) 

z 2 +n 2 dz^dz' 


(7.60) 


In applying this to Equation 7.59 some cancellation occurs leaving 


12 2ttz 


Z cos(n^) 
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+ I n (naz) S(z o ) 
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dz. 


(7.61) 

The reason for obtaining this form is to permit integration by parts, the result of 
which is 
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(7.62) 


The first line comes from the integrated terms evaluated at the field point z; the 
second line comes from these terms evaluated at the extremes of the integration 
range; and the last line contains the remaining integrals. 

The first line can be simplified by noting that the expression in square brack 
ets is the Wronskian W[K n (x) , I„(x) ] - 1/x. In anticipation of the final result, 

this term will be called a^ w so that 


where 


a 12 Q 0 U n0 


Q n W ” TT~ s ( z ) I cos(mB^) 
u z,rz m=l 


(7.63) 


(7.64) 


and 


a” =naz 2 S(z 2 ) I n (naz)K^ (naz 2 ) - naz x S(z x ) K^naz)!^ (nazj) 


- K^naz) naz Q ^~ I„ (naz o )dz Q - I n (naz)| naz Q dz o ° ^ (naz ° ) dz ° (7-65) 


The more interesting term here is c^ w with its non-decaying cosine series. 
Because of its importance, two independent ways of interpreting this term will be 
given. The simpler scheme is to rewrite the series as follows. 

Y cos(mB^) = i X cos(mB^) - i (7.66) 

m=l z 


The double sided series is related to the delta function series 


£ cos(m$) — 2n S ($>) (2.67) 

m=-® 

Since this is 0 everywhere except at $ = 0 , it is also 0 in the limit as 
approaches 0. Similar analysis of the series in Equation 7.66 shows that it too is 

0 except for £ - 0, 2 tt/B, 4tt/B, etc., which is to say on the trailing vortex, 
sheets. Since it was emphasized above that the downwash formulas are valid in the 
limit as the field point approaches the vortex sheets, we conclude that the correct 
value for the sum in Equations 7.64 and 7.66 is -1/2. Thus, 
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( 7 . 68 ) 
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The conclusion of this analysis of Equation 7.64 is that an expression with 
only oscillating terms and no mean value is correctly interpreted as having a mean 
value and no oscillation! In fact, this term represents the mean (circumferentially 
averaged) swirl and axial induction in the wake that seemed to be missing and thus 
resolves the paradox described in conjunction with Equations 7.48 and 7.51. The 
significance of the minus sign in Equation 7.68 is that the downwash for positive 
blade loading is, in fact, down. Comparison of this expression with classical 
propeller momentum theory will be given In Section 8 and the alternative means of 
deriving Equation 7.68 will be given in the next subsection. However, first we take 
time to comment on the results derived so far. 

It was stated above that Equations 7.43-45 constitute the final form of the 
propeller kernel. This is indeed the form used for the discretization described in 
Section 9. However, both the wake and bound terms of the kernel are singular where 
the load point radius z o and field point radius z are equal. The leading singu- 
larity in the wake term has just been discussed and It is clear that Equation 7.44 
is not suitable for numerical work for small z-z Q . The scheme used in the com- 
puter code is to divide the radial integration range into 3 sections: a narrow 
singular section (typically Az o = .02) centered on z, an inboard section extend- 
ing from the root to the singular section, and an outboard section extending from 
the singular section to the tip. In the singular section a special form is used for 
the wake portion of the kernel based on the integration by parts scheme presented 
above. In fact, for best results, another Integration by parts is used to deal with 
a logarithmic singularity driven by the second derivative of S(z o ). 

The bound term of the kernel has a logarithmic singularity at z o = z that is 
also manifested in the non- convergence of the harmonic series. This too is dealt 
with by dividing the radial integration range into three sections. For the singular 
section, a patch is made to a special form of the wing kernel. The logarithmic 
singularity is then integrated by analytical curve fitting. The special form of the 
wing kernel was derived by the same pressure potential method used for the propel- 
ler. The resultant form is similar to Equation 7.44 with the Fourier transform but 
without the summation. 

Integral of Body Force Term in Momentum Equation 


We now return to the singular behavior of the wake term and approach it by a 
different method. This is a variation of the scheme used by Lordi and Homicz (Ref- 
erence 7) in their pressure potential analysis of ducted rotor flow. They studied 
the singularity using the velocity potential rather than the downwash expression, 
but the concept is similar. The idea is to include the body force term in Integra- 
tion of the momentum equation. When this is done, it is found that the singulari- 
ties of the force term and those of the pressure term cancel exactly, as will be 
seen. 

Representation of the body force for steady loading was given by Equation 3.20 
(with w Q = 0) in coordinates fixed in the fluid. The lift component is 

f L - - 5 (£+FA) f L (y+Ut) (7.69) 
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f is the force per unit volume on the fluid whereas f L is the force per unit area on 
the helicoidal surface. 7 is the coordinate running backward along the helical 
advance path as shown in Figure 4 so that 7 +Ut represents the blade section motion. 
For convenience in manipulation, Equation 7.69 can be expressed as a delta function 
integral as follows . 

f L = - 5 (£+FA) Jf L (7 0 )S(7-7 0 +Ut)d 7o (7.70) 

By using the same changes of variable that were described in conjunction with Equa- 
tions 3.13 to 3.15, the volume force can be re-expressed as 

f L - - y 5 (^-^-^ FA )/fL(7 0 )5(-^x-^FA-7 o +Ut)d7 o (7.71) 


When this is changed to fuselage -fixed coordinates via x— x^+Vt and then the sign 
convention for x is changed to positive downstream via x^-x, time drops out, 
leaving 

f L - - § S(hJf L (7 0 )S(™-^FA-7 o )d7o (7 ' 

where % was given by Equation 7.31. By expanding the first delta function in a^ 
series and applying the property 5 (ax)- 6 (x)/a to the second delta function, this 
can be further modified to read 

f L - - 3?r 1 e imB ^ Jf L ( 7 Q )S(Ax)dr o 

■4 ' 1 h f^=_C0 

The number of blades effect has been accounted for with the appearance of 
a coefficient and in n=mB as in Section 3. In coefficient form, the lift 
unit area 

£j(\) - 

If we call the associated force/unit volume f L , then the acceleration due 
body force becomes 

f f L - iSr l e imB3 ;AC p (7 o )5(Ax)dr o (7.75) 

H Q H71L m =-oo 

This was the term avoided in Equation 7 . 3 by integrating above or below the airfoil 
as suggested by Figure 10. But, Equation 7.75 can now easily be Integrated using 
the method of following the particle outlined at the beginning of this section. 

Since f L is the force normal to the helicoidal lifting surface, its integral is 
the disturbance velocity normal to the surface, w. Expressed as an angle, the 
result is 


(7.73) 

the B as 
force per 

(7.74) 
to the 


g - a = gi J a e imB ^aAC p ( 7 o )H(Ax)dr o (7.76) 

Since the integral in this equation Is the source S defined by Equation 7.58, the 
downwash due to the body force can be rewritten as 

a = S(z) H 1 + 2 £ cos(mB^) (7.77) 
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Now it can be seen that, when the downwash caused by the body force is added to that 
caused by the pressure term in the momentum equation, the series in Equation 7.64 is 
cancelled term by term for m=l to «. Only the zeroth term in Equation 7.77 remains 
so that the circumferential average term is -BS(z)/47tz, as was found above by the 
limiting process. 

Both the limiting process and the body force cancellation schemes are valid 
means for dealing with the singularity at z o - z and finding the circumferential 
average induction. Because of the subtleties involved, it is worthwhile to under- 
stand both methods. Furthermore, it will be seen in the next section that the body 
force method is needed to understand the circumferentially average axial induction 
and compare it to classical momentum theory. 
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SECTION 8 

RELATIONSHIP BETWEEN LIFTING SURFACE THEORY 

AND 

CLASSICAL PROPELLER MOMENTUM THEORY 


In the previous section, equations were developed for the downwash normal to 
the helicoidal advance surface caused by a specified distribution of lift pressure 
on the surface. This was in complete analogy to the wing downwash equations whose 
advance surface is a plane aligned with the flight direction. The equations are in 
the form of a Fourier series whose zeroth term gives the circumferential average of 
the downwash (or induction) and whose higher terms give the blade- to -blade varia- 
tion. In this section, we concentrate on the circumferential average terms. They 
will be formulated to give axial and tangential induced flows as functions of ele- 
mental thrust and torque. This is as opposed to the treatment In previous sections 
that gave downwash as functions of lift force. The purpose of this is to compare 
the lifting surface theory with classical propeller momentum theory, which deals 
only with circumferential average effects. In the process, insight is given as to 
the roles of the various terms, the effect of compressibility , and a method for 
treating the non-linearity associated with finite axial induction. The derivation 
is given first for the general case of a rotor with loading distributed along the 
chord, i.e. a "thick" rotor. The blades can have sweep and offset. Then, for 
comparison with momentum theory, the loading will be specialized to an actuator disc 
form. 


Mean Induction for "Thick" Rotor 

The derivation for the axial and tangential induction components is based on 
the integration of the linearized momentum equation 


Dy 

Dt 


i-(F-Vp) 

O 


( 8 . 1 ) 


The integration procedure is similar to that in Section 7, but different in detail. 
We use the result, established at the end of Section 7, that the zeroth term of the 
pressure field can be taken from Equation 7.8 and the zeroth term of the body force 
can be taken from the original source representation. We start first with the 
pressure field. 

The circumferential average of the pressure is given by the n=0 term of Equa- 
tion 7.8 

p ^ | U ° f ^ ^ ^F( k x) J o( aZ < K ) H 0 1) ( aZ > K ) dwdZ o < 8 ' 2) 

Z h “® 

The radial wavenumber K is pure imaginary for n=0 and q=0 

K = = f/jV (8.3) 

so that the Bessel function product becomes 
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(8.4) 


JH = y ~ Iq ( a ^ Z< | w | )K 0 (a/Jz., | w | ) 


where I 0 and K 0 are modified Bessel functions, as shown in Section 3. 
is the Fourier transform representing the axial and tangential force distributions 
as given by Equation 3.69. It is important to note that, for n=0, there is no 
contribution from the tangential force to the circumferential average of the pres- 
sure. As will be seen, the consequence of this is that there is no mean swirl 
upstream of the rotor. To return to the original force representation from the 
Fourier transform, Equations 3.69, 3.61, and 3.21 are applied in that order with the 
result 


v**> e 





(7 )e 


dr„ 


(8.5) 


where f x is the distribution over a blade of the axial force per unit blade 
area. As a Fourier transform, the integration range runs from -« to °o but, of 
course, the blade loading is non- zero only between the leading edge and trailing 
edge. Substitution of this expression into Equation 8.2 gives the pressure field 


p - }}f x (7 o )|^I 0 (a^z < |w|)K 0 (a^| w |)e- ia ^ dwdr o dz o (8.6) 

where Ax, the axial distance between source point and field point, is repeated here 
from Section 7 


Ax - - -2 - 

r T a o 


2 az 


o FA 
D 


The gradient of p involves derivatives in the x, <f> t and z directions, 
see that the gradient in the tangential direction is zero and, in the present 


(8.7) 

But we can 


linearization, we will ignore the radial gradient. The remaining component, idp/dx, 
is trivial to find. The operations to compute the axial velocity disturbance asso- 
ciated with the pressure (or bound) term are symbolized by 


v" 


-f * 
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3p(x') dx' 

ax' 


v 


( 8 . 8 ) 


The V under the integral requires some discussion at this point. It represents 
the velocity of the particle whose acceleration is being integrated as a function of 
particle position x' . Equation 8.8 could be integrated numerically by setting 
V=V o +v xo and updating V at each step in x' . However, to deal with Equation 8.8 
analytically, V is set equal to its "average” for the integration range, namely 
V(z), and brought outside the integral. Interpretation of the word "average" will 
be given below in the momentum theory discussion. 


For the integration in Equation 8 . 8 , there is no contribution at w=0, as there 
was for the higher order terms in Section 7, because the limit of the integrand is 0 
at w=0. The impact of differentiating and then integrating with respect to x is 
simply to change the constant in Equation 8.6 with the following result 
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(8.9) 


for the velocity disturbance associated the pressure field. From consideration of 
symmetry in the u integral, this can be written in pure real form as follows 
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To find the velocity disturbance associated with the body force, we pick out 
the zeroth term from Equation 7.73 and interpret it for both the axial and tangen 
tial force components as follows 


f x = 2fcf f x(7 0 )S(Ax)dr 0 

and 

f 4 - 2fcr4<v 5 < ax > dr ° 


( 8 . 11 ) 


( 8 . 12 ) 


The minus sign that was in Equation 7.73 has been dropped because of the following 
convention. f x is the axial force per unit volume acting on the fluid and f x is the 
axial force per unit area acting on the blade surface. The usual sense of positive 
thrust for a propeller is associated with a downstream (positive x direction) force 
on the fluid; thus, the minus sign is not needed here. Similar remarks apply to the 
tangential force. The integral to find the velocity caused by the body force is 


This is easily applied to Equation 
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8.11 with the result 
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(8.14) 


where the V has been brought outside the integral as V for the same reasons 
presented after Equation 8.8. The same procedure applied to Equation 8.12 gives 
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(8.15) 


This completes the derivation for the circumferentially averaged induction for 
steady loading. Tangential velocity is given by Equation 8.15. This simply states 
that v^ 0 is proportional to the tangential force integrated up to the axial field 

point x. The axial induction is the sum of the bound and wake terms 
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where, as with the tangential component, the wake term is proportional to the inte- 
grated force up to the field point. The bound (or pressure) term in Equation 8.10 
is in the form of a Fourier integral and contains the compressibility effect via the 
P in the arguments of the modified Bessel functions. The w integral could, per- 
haps, be evaluated in terms of elliptic integrals, but this doesn't appear to be 
particularly useful. A better approach is to move the T 0 integration to the inside 

and compute the Fourier integral of f x . In the discretization process described in 
Section 9, the loading is expressed as a series of shape function with unknown 
coefficients. The shapes are chosen such that their transforms can be computed 
analytically. This leaves only the w and z o integrals, which can be computed 
rapidly. 


Actuator Disc Results 

The propeller loading can be collapsed onto a zero thickness disc by setting 
the face alignment FA to 0 and replacing the force-per-unit-blade -area representa- 
tion by thrust- and torque-per -unit-radius as follows 


- fj |r «<V 

- rh f- t( v 


(8.17) 

(8.18) 


Dealing first with the tangential velocity, we substitute Equation 8.18 for the load 
representation into Equation 8.15 and perform the r o integration. The result is 
simply 


V ' SS Vr* f H(ax) (8 - 19 > 

where here the - signifies actuator disc results. This shows, via the step func- 
tion H(Ax) , that tangential velocity component jumps from zero ahead of the rotor 
to its final value immediately behind the rotor. The behavior is shown in the top 
sketch of Figure 11 and agrees with textbook treatment of propeller aerodynamics 
(ref. 31). 

For the axial induction, we substitute Equation 8.17 into Equation 8.10 and 
perform the F o integration with the following result for the bound or pressure 
term 


V xo 


Y^'p\r 2 T f ^ I wI oC^)K 0 (^z > w)sin(w -^) dw dz 0 


( 8 . 20 ) 


And from Equation 8 . 14 the wake term is 


Vxo 2-nzp Vr| dz 


€ H(Ax) 


( 8 . 21 ) 
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The total axial induction for the actuator disc approximation is the sum of these 2 
terms 


The wake term of the axial induction, shown in Figure lib, jumps at the rotor disc 
in the same fashion as the tangential velocity. The bound term is odd in x and has 
the general behavior sketched in Figure 11c. 

Since the bound term of the axial component also jumps at the disc plane, we 
will evaluate the magnitude of the jump to see how the 2 terms work in concert. 
Because the argument of the sine function in Equation 8.19 includes the product wx, 
it is clear that the behavior near x = 0 is controlled by the behavior of the inte- 
grand at large w. Thus, we can use the large argument asymptotic form for the 
Bessel functions to evaluate the inside integral in Equation 8,20. 


I - J'“wI 0 (^wz < )K 0 (/3wz > )sin(w^)dw 

Upon substitution of asymptotic forms from Reference 25, this becomes 

X = \ f e z cJ sin(w ^-)dw 

2fS^zT 0 J o r x 


(8.23) 


(8.24) 


which is given by a standard integral as follows 


I 




(8.25) 


Now, in the remaining z Q integration, any section of the integration range not 
containing the field point z will produce 0 in the limit as x-*0. However, we still 
must consider a small band centered on z to find any contribution from that point. 
Thus , we evaluate 


vt(0 + ) 


dT 


Lim 


hn z p V^r|z dz x-o 


x 

r T 


.z+Az 


z-Az 


dz„ 


(^) +0 Z ( Z - Z o) 2 


(8.26) 


The term in square brackets is 


[] = i tan- 1 ^ zj| 


Az 

-Az 


(8.27) 


If Az is held fixed and finite, the limit of this is n/f) and the result for 
Equation 8.26 is 
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(8.28) 


v" o (0 + ) = 


-1 


dT 


4ttz p Vy9 z r z dz 


and the jump in the bound term at x=0 is twice that or 


Av 

xo 


-1 


dT 

2?rz p Vf3 Z Tj ^ z 


(8.29) 


The jumps in axial velocity from the bound term in Equation 8.29 and the wake 
term in Equation 8.21 can now be compared. For incompressible flow (y3=l) , the 
jumps exactly cancel so that the axial flow is smooth through the disc as in the 
classical propeller momentum theory (ref. 31) and sketched at the bottom in Figure 
lid. However, for compressible flow, the jump in the bound term is larger so that 
the sum of the 2 terms is discontinuous as also shown in Figure lid. The velocity 
jump is required for mass continuity because the density in front of the rotor is 
lower than that behind. This was recognized by Vogeley (ref. 32) but could not be 
included in his actuator disc analysis. 


Sample Calculations 


The discontinuous jumps in velocity found in the preceding section are, of 
course, artifacts of the actuator disc representation of the loading. With the more 
general "thick rotor" formulation above, the velocity variations are smooth but 
still influenced by compressibility. Here we present some sample calculations of 
axial induction for realistic rotor geometry and loading using Equations 8.10 and 


8.14. 


The study rotor is representative of the front rotor of a counter rotation 
Prop-Fan designed for a flight Mach number of 0.72 and an advance ratio of 2.94. 
Figure 12 at the bottom gives radial distributions of chord- to -diameter ratio and 
thrust coefficient. At each radial station the chordwise distribution of loading is 
assumed to be constant. The top of the figure shows the axial induction profiles 
at several axial stations. Note that behind the rotor the induction is forward 
outside the blade tips. This obviously could have a strong effect on a rear rotor 
if the streamline from the front rotor contracts inside the rear rotor tips. A 
trailing wing would also be affected. 

Figure 13 plots the axial induced velocity normalized by the flight speed at 
radius ratio z=0,78 as a function of axial coordinate. Incompressible calculations 
for the same case are also shown. As described above, the effect of compressibility 
is to increase the induced velocity in front of the rotor and reduce it behind. 

There is little effect at blade mid-chord. The compressible actuator disc results 
above gave a discontinuous axial velocity whereas the thick rotor theory predicts a 
smooth distribution. 

The compressibility effect just shown has an interesting influence on chordwise 
blade loading. Higher Induction in front tends to unload the leading edge and lower 
induction behind tends to add loading to the trailing edge. This result could not 
be deduced from lifting line theory, which deals only with mid-chord induction 
values. The consequence is that compressibility produces the need for more camber 
in the blade sections. The effect is in addition to the aft movement of the loading 
on 2D airfoils caused by non-linear compressibility. 
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Comparison with Prope ller Momentum Theory 
and Effect of Non-linearity 


Simple propeller theories require the axial induction at the blade midchord so 
that 2D airfoil data can be used. We will call this midchord value v xd and note 
from the figures just discussed that it is 1/2 the ultimate wake value. In Weick's 
book (Reference 33), Newton's second law takes the following form. 


dT = ( dA * pV ) * ( 2v xd ) 

(mass/unit time) * (velocity increase) 


(8.30) 


where dA is an annular element of disc area and V is the velocity through the 
disc. Thus, the velocity at the disc can be expressed 


v. 


xd 


1 dT 
2p V dA 

O 


(8.31) 


To compare this with the theory of the present report, note that Equation 8.14 shows 
the wake velocity to be 


2v. 


dT 


xd 


2nzp VrS dz 


(8.32) 


When it is recognized that the elemental area dA = 25rzr^dz, then it can be seen that 
the circumferentially averaged induction of the present theory is exactly the same 
as classical momentum theory provided that V has the same meaning. This is 
discussed in the next paragraph. 


In Equation 8.30, V is related to the mass flux through the disc and must 
therefore include the flight speed plus the axial induction at the disc v xd . Thus, 
we can write 


1 dT 

V * d “ 2p o (V o +v xd ) dA 


(8.33) 


This is the form of the induction equation shown in many textbooks. Since v xd and 
thrust are not known a priori, loading is usually found by iteration. For complete- 
ness, note that these same procedures applied to Equation 8.19 lead to the standard 
relation between torque and tangential induced velocity 


v - 1 ^ 

2p Q (V 0 +v xd )r dA 


(8.34) 


In the original development of the author's helicoidal lifting surface theory 
(ref. 15), the derivation paralleled that for wing pressure potential theory. In 
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the wing problem, there is no axial induction and it is correct to set V equal to 
the flight speed V Q . This same linearization works well for propellers in high 
speed flight; but, at takeoff conditions, V o and v xd can be about the same and the 
linearization gives poor results. For the static condition, it fails completely. 
Then Equation 8.33 must be iteratively solved with v xd included on the right 
hand side of the equation. 

Now that a means to deal with the axial induction non-linearity has been estab- 
lished for the circumferentially averaged flow, it still must be decided how to deal 
with that effect in the full lifting surface theory. The approach taken is that we 
force the non-linear momentum equation (Equation 8.31) be satisfied in the circum- 
ferential average sense at each radial station by iteration. This establishes the 
value of V as a function of z. Blade -to -blade variations in induction and 
unsteady effects are calculated using equations in Section 7 with the flight speed 
set to V. Effectively, the solution is then linearized about the non-linear 
momentum value at each radius . The iteration scheme is discussed at the end of 
Section 11. 
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SECTION 9 

DISCRETIZATION OF INTEGRAL EQUATION 


This section develops the method for computing blade loading from the downwash 
theory derived in Section 7. Blade loads are needed as input to the noise equations 
in Sections 4 and 5 and to the performance equations in Section 11. The downwash 
equation was given by Equation 7.33, which is repeated below. 

« - /;aC p (7 0 > z 0 ) K dr o dz o (9.1) 

This gives the downwash angle a as an integral over the blade surface of the coef- 
ficient of lift pressure times a kernel function. The kernel function is given by 
Equations 7.43 and 7.44. The difficulty with Equation 9.1 is, of course, that the 
unknown pressure is under the integral whereas a is known on the blade from 
geometry. Specifically, a is the angle the flow must be turned by the loading so 
as to be tangent to the camber surface. To solve Equation 9.1, it is inverted by 
discretizing the lift pressure (i.e. writing it in the form of a series of load 
elements with coefficients to be found), converting it to a matrix equation, and 
inverting the matrix. The procedure is similar to that described in Reference 16. 
The main difference is in the refinement of the loading elements from a quasi -vortex 
lattice representation to true lifting surface. 

In the next paragraph, the general procedure is outlined. Then, the inversion 
scheme is worked out in detail. 


General Form of Inversion Method 

First, we write Equation 9.1 in generalized form as 

a = f AC K dA (9.2) 

j p 

where dA is an element of blade area. Flow tangency will be enforced at a number of 
control points on the blade surface counted by the index /*. The array of downwash 
angles at the control points forms the downwash vector so that Equation 9.2 

becomes 


- J AC P K „ (9 ' 3) 

where the kernel function is now evaluated at the control points. Symbolically, the 
discrete form for the loading can be written 

AC p - l L v S„ (9.4) 

V 

where the S^' s are a family of shape functions of unit amplitude and the L u 's 

are their coefficients which are to be found. Inserting this loading representation 

into Equation 9.3 gives 


W„ = ; l K s v dA (9.5) 

V 

Now the integral is moved inside the sum with the result 

W M - l L„ J S„ K„ dA (9.6) 

V 
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The integral is a function of both indices and becomes the matrix element 

V - ; S„ dA (9.7) 

Thus, by discretization, the integral equation has been converted to the matrix 
equation 


W - I L, K (9.8) 

U 

Its inverse is 

K “ l (9.9) 

ft 

With this system, the known downwash can be used to compute the coefficients L u 
and these can be used to construct the loading with Equation 9.4. 

This has been a schematic description of the inversion process. The actual 
method described below for the propeller equation is more complex, but follows the 
same principles. 


Inversion of Propeller Integral Equation 

The field point in the kernel given by Equations 7.43 and 7.44 is given in 
cylindrical coordinates. In this section, the control point locations must be 

expressed in terms of blade surface coordinates so we switch x and $ to heli- 
coidal surface coordinates using the Inverse of Equations 3.11 and 3.12, namely 

x = u 7 ' tP = h - < 9 - 10 > 

r *i - - if 7 - & = - U ( 9 - n > 

Note that the sign convention for x has been changed from that in the acoustics 
sections to positive downstream for the aerodynamic applications starting with 
Section 7. With this transformation, the angle between the source and control point 
helicoids becomes 


a i °o FA 
2 ^ r T 

and the axial separation between source and control points becomes 


(9.12) 



+ Ax e 


where Ax^ is shorthand notation for 


Ax 




az £ az o FA 
a rj " a o r T 


(9.14) 


(9.15) 


With these changes In notation, the appropriate expression for downwash is obtained 
by combining Equations 7.42, 7.43, 7.44, and 9.1 
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x I[K f0 K e JH^w) - K f0 K ? 0^(0) ] dw dz c 


(9.16) 


The system for load paneling and control point locations is shown in Figure 
14. It shows the blade divided into NCP spanwise panels which are tapered at con- 
stant percent chord. Control points, where flow tangency is enforced, are in chord- 
wise arrays at radii z. . The index system for the control points is the same as 

in Reference 16, namely that i-1 , 2 NSM counts the radial locations and m 

_ ]_ i 2 NCP counts the chordwise positions. Downwash velocities are given by 

the vector whose index combines the chordwise and spanwise counters as follows 


/i = (i-l)NCP + m 


(9.17) 


Downwash is not computed simply at the control points but rather is averaged over a 
small chordwise range by using the following operation 

w„ = J a G(r-r i ,z i ) dT < 9 - 18) 

G is a weighting function that centers the averaging at the chordwise location T- at 
radius . This chordwise location is given by 


r- = 


MCA 


-B d + (m - 1 +CP)A 0 


(9.19) 


where MCA is the sweep defined by Figure 4, B D is chord to diameter ratio, CP is 
the location of the control point within the panel (0.5 for the middle), and the 
panel width is 


2B d 

NCP 


Various forms for G could be used but the best seems to be 


G(r - r s- z i ) exp 


c 


zA^ ( r-r-) z ] 


(9.20) 


(9.21) 


This is a Gaussian function with unit area. Dependence on z i comes in via A, the 
width of the Gaussian function at its half amplitude points. This is tied to the 
local panel width on the blade by the formula 


2 B 

— 2. * CHW = A * CHW 
NCP 


(9.22) 


where B n is chord- to -diameter ratio and CHW stands for control point 1/2 width. 

A is the width of the weighting function at its 1/2 amplitude points normalized by 
the width of a load panel. Averaging the downwash in this manner is reasonable from 
a physical point of view; however, the main reason for doing it, as will be seen, is 
that it damps the w integrand at large values of w and guarantees convergence of 
the integral. Furthermore, on blade sections with supersonic speed, the kernel 
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function has very strong gradients near Mach waves. When control points are near 
these Mach waves, results can be sensitive to the exact location of the control 
points if the kernel is not averaged. 

The scheme for load elements is also shown in Figure 14. As mentioned above 
the blade is divided into NCP spanwise panels of constant percent chord. Spanwise 
variation of loading within panels is given by a series of shape functions R as 
sketched at the upper right in the figure. Each of these has built into it an 
elliptic falloff at the tip since this is the correct behavior for wing tips. A few 
different families have been evaluated but the one currently used by default in the 
computer program is 


R j “ z o sin [j c °s" 1 (z 0 )] (9.23) 

The first 7 of these functions are plotted in Figure 15. The family was adapted 
from one commonly used in wing methods by multiplying by z Q to drive the functions 
to 0 at z Q - 0. It can be seen that the functions vary most rapidly at the tip 
where pressure gradients are expected to be the greatest. 

Chordwise variation of the load elements is also shown in Figure 14. A series 
of overlapping elements C s is used with some special treatment at the leading and 
trailing edges as shown at the lower right in the figure. At the leading edge, the 
user of the computer program can choose the triangular element shown with the dashes 
or the preferred singular element shown by the solid curve. The latter element is 
formulated so that, with the adjacent triangular element, it gives the correct 
leading edge behavior for 2D cases. At the trailing edge, the basic element goes to 
for subsonic trailing edges [M rel cos(trailing edge sweep angle) < 1] satisfying 
the Kutta condition. For supersonic edges, a finite jump is allowed since the Kutta 
condition doesn't apply. There is a blending at the sonic radius between the 2 
types of element. Chordwise elements C 5 are counted by h which runs from 1 to NCP. 
The combined index for the spanwise and chordwise variation is 

u = (j-l)NCP + n (9.24) 

The associated loading can now be written 


ac p ( V z o> “<r Ik R j< z o) c 5 (r o ,z 0 ) 


J,n 


(9.25) 


Note that this distribution is continuous in the chordwise direction and continuous 
and smooth in the spanwise direction. This degree of smoothness is needed for 
supersonic flow because any artificial discontinuities in loading can cause Mach 
waves that lead to unstable solutions. If the method were only to be used for 
subsonic flow, it would be much less sensitive so that a crude load system could be 
used, perhaps even a system of constant pressure panels. 

The load elements C s are each different in the sense that they appear at vari- 
ous values of T 0 , depending on blade geometry. By using a change of variable, they 
can be expressed in terms of 3 common functions as follows 


c s( r o. z o) = F 4 <r 0 -r .) 


(9.26) 
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(9.27) 


r - is the helicoidal distance to the leading edge of element h 

on 


r 


on 


MCA 


B d + (m - 1)A 0 


Figure 16 defines the 3 types of load element in terms of the panel width A 0 . 

We can now insert the discretized loading given by Equation 9.25 into the 
downwash equation (9.16) and perform the control point averaging indicated by Equa 
tion 9.18. The result is the desired matrix equation 


- I K V 

! v 


(9.28) 


where the matrix elements are, using the crossed integral sign defined below, 


K 


‘fiU 



and the chordwise- integrated kernel is 


(9.29) 
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(9.30) 


In Equation 9.29 the significance of the crossed integral is that integration of 
z through the control point radius z t does not exist in the ordinary sense. It 
must be interpreted as the "finite part of an infinite integral" or the Mangier 
principal value 27 because of the singularities discussed in Section 8. 

If we now shift the integration variables with r-+r-r s and r o -»r o - r o - , then the 
chordwise - integrated kernel becomes 


K~ 
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tqa( ^ ^ + Ax-) 


I- l e^ K 50 K f I n ( |n+q | az < )K n ( | n+q | a Z> ) 


0111 47T(7 2 ZZ„ 
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e wt K fO K f 


- K- 0 K-JH n (0) ] dw(9 . 31) 


where the double chordwise integrals are 
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Ik = /G-(r)/F A (r o )H(Ax')e iaq(? " Q> dr o dr 


(9.32) 


and 


- — ) 

n s - J G s (r)/F A (r 0 )e ° dr o dr 


The argument of the step function is 

- (tt= • “ ■ [ 
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FA 

° a o r T J 


(9.33) 


(9.34) 


Finally, we can define the axial displacement between the control point and the load 
element as 






(9.35) 


and arrive at the final form for the chordwise- integrated kernel 
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e ^(w) - 


K C0 Kf 


0^(0) ]dw (9.36) 


The first and second terms give the wake effect and the bound effect, respectively. 

Major results of this derivation can be summarized as follows. The matrix form 
of the lifting surface integral equation is given by Equation 9.28. The kernel 
matrix is to be computed from Equations 9.29 and 9.36, The inverse equation (9.9) 
can be solved for the loading vector L u given the vector of downwash angles 
W . The distribution of lift pressure over the blade surface can be constructed 

from Equations 9.25 and 9.26. Details for computing the downwash vector are given 
in Section 10. 


Comments on Evaluation of Kernel Function 


To evaluate I--, note that the integration ranges of r and T are very small 
because of the local nature of the G and F functions. Thus, the exponential can be 
set to 1 for q "not too large". This amounts to neglecting phase variations within 
one load element; phase variations from panel to panel are still included via the 
exponential containing x . If this approximation is not satisfactory for a given 
set of panels, then the number of panels can be increased to make the G and F func- 
tions more local. Thus, for the integral in the wake term, we have 

- jG s (D/F A (r o )H(Ax') dr o dr (9.37) 

which is easily evaluated by numerical methods. It exhibits the wake behavior 
explicitly since it varies smoothly from 0 upstream of the rotor to 1 behind. 
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For evaluation of II- , note that Equation 9.33 can be written as the product of 
2 integrals 

n X * X, (9.38) 

mn i 


The first, 

x = jG-(De' ialu ' 5) ' dr (9.39) 
is the result of averaging over the control point range and can be evaluated analty- 
ically with the result 

(9 ' 40) 

where A is defined by Equation 9.22. As mentioned above, this has the effect of 
damping the integrand in Equation 9.36 and guarantees convergence of the integral. 
The second of the 2 integrals in Equation 9.38 is 




JX< r o) 


ia(v-q ) 
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dr o 


(9.41) 


This is the Fourier integral of the load element shapes in Figure 16 , which can also 
be evaluated analytically except for the singular leading edge element. For 
example, for the mid- chord elements 
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(9.42) 


This is pure real because F A is an even function. The leading edge element must 
be integrated numerically. This is done once during a computer run and stored as an 
array for later interpolation. The cost is negligable. 

In Equation 9.36, the integral over w is also a Fourier transform. Since it 
is needed at several values of the variable x^, it is best evaluated by fast Fourier 

transform (FFT) methods. The integrand is conditioned for good behavior at the 
origin and the infinite limits cause no problem because of the damping provided by 

the X function. 

Finally, some remarks are needed on the radial integration denoted by Equation 
9.29. As was mentioned above, the double bars through the integral sign signify 
special treatment as the load point passes through the control point, i.e. near 
z Q = Zi . The wake term and the bound term are integrated separately to keep 
control over their singular behavior. For the wake term, the radial integration is 
divided into 3 ranges: a narrow range centered on the singular point and 2 outer 
(non- singular) ranges covering the remainder of the blade. The non-singular ranges 
are done with a direct integration of the first term of Equation 9.36. In the 
singular range, the integrand is reformulated using integration by parts as 
described in Section 7. This exposes singularities of order (z i -z Q ) , (z-zj , and 
log| Zi -z 0 | which can be integrated analytically because the mode functions and 
their derivatives are given in analytic form. 

Radial integration of the bound term in Equation 9.36 is also divided into a 
singular range and 2 non-singular ranges. In the non-singular ranges the second 
term is integrated directly. For the singular range, the order of the singularity 


- 73 - 


had to be established. This was first accomplished by some very tedious series 
expansions of the integrand using asymptotic forms of Bessel functions to expose a 
logarithmic singulartiy. The result was confirmed by analysis of the bound term of 
the steady wing kernel in Equation 7.39. It was found that, when this is integrated 
in the chordwise direction first, as we have for the propeller kernel, the remaining 
spanwise behavior is logarithmic at chordwise locations within the load element. In 
fact, strength of the singularity is proportional to the chordwise derivative of the 
load function. Since the behaviors of wing kernels and propeller kernels are 
expected to be the same near the singularities, it was decided to patch to a wing 
formula in the singular region. Unfortunately, Equation 7.39 applies only to steady 
loading and the published unsteady formulas are unwieldy. A suitable form of wing 
theory was derived by applying the methods of Sections 3 and 7 to translating 
sources, including the procedure of separating the wake and bound terms. The result 
is again a Fourier integral involving a Bessel function but without the infinite sum 
that appears in the rotating case. It was verified by numerical experimentation 
that this form of the bound wing kernel matches the propeller kernel close to the 
singularity. Because of its simple form, the bound wing kernel can be computed with 
good enough accuracy to integrate it by numerical curve fitting. This is the scheme 
used in the computer program. 
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SECTION 10 

BOUNDARY CONDITIONS AND INVERSION OF INTEGRAL EQUATION 


In the preceding section the lifting surface integral equation that was derived 
in Section 7 was discretized for conversion to a matrix equation as follows. 






'\iV 


( 10 . 1 ) 


is the vector of downwash angles on the blade surface at control points which 
are counted, as shown in Figure 14, by the index is the vector of loading 

coefficients whose index u corresponds to members of a family of loading element 
shape functions. K is the matrix of influence coefficients giving the downwash at 
control point /i due to a unit load from element v. Equation 10.1 is to be solved 
by matrix inversion 


I K 


1 W 


( 10 . 2 ) 


giving the loading vector for specified downwash angles at the control points. 

Note that Equation 10.1 gives the flow deflection due only to the blade load- 
ing. There are other deflections, due to blade thickness and nacelle blockage, that 
must be accounted for before the matrix equation can be applied. The "air angle 
accounting" system is the subject of this report section. 


Boundary Conditions 

Definitions of the angles to be used are given in Figure 17. The basic 
requirement for flow tangency is that the air angles must match the hardware angles 
at the control points. As shown in Figure 14, control points are placed at radius 
ratios z A , where the radius index i runs from 1 to NSM, and are placed one per 
panel across the chord with the index m. The combined index for matrix rows is 

H = (i-1) * NCP + m (10.3) 

Hardware angles are shown in the figure at the upper left. These are specified 
in terms of a blade angle (or twist angle) 9 h and camber angle 0 cam . As described in 
conjunction with Figure 4, the angles are measured in the cylindrical surface, z i 
= constant. Thus, to construct the hardware angles the blade is "cut" by the cylin- 
drical surfaces z i and rolled out flat before the angles are determined. Blade 
angle is a function of index i only. It is determined by constructing the chord 
line and measuring its angle from the plane of rotation. Camber angle is a function 
of both the i and m indices. It is computed according to the sketch at the 
bottom of Figure 17 from 


„ „ dh 

P cam dX 

at chordwise positions 

y m - 1 + CP 
S = NCP 


(10.4) 


(10.5) 
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where CP is the non-dimensional position of the control point within the panel. The 
panel midpoint corresponds to CP=0 . 5 . In Equation 10.4, the small disturbance 
assumption was invoked to approximate the tangent of the camber angle by the angle 
itself. 


Flow angles are defined by the velocity triangles in Figure 17. Forward flight 
speed is V o and the rotational speed at radius ratio z i is 7rNDz i . Thus, the 
section advance angle at radius z i is 


tan 6 

o 


?rNDz i 


advance ratio 



( 10 . 6 ) 


For the blockage, a separate calculation is performed with any appropriate axisym- 
metric flow program to represent the effect of the nacelle, hub, and spinner. Axial 
velocity components are stored on a grid that encompasses the propeller position and 
then interpolated to the blade control points to find the blockage angle from 




block 




(10.7) 


as a function of both indices i and m. Finally, the thickness effect Is 
computed from the theory in Section 3 with details given at the end of this section. 


Using the notation of Figure 17, the match of the hardware angles to the flow 
angles can be expressed for each control point by 


a. + a-. + <j> + <f), 


block 


7t> + 0 

n cam 


( 10 . 8 ) 


or, in terms of the unknown angle, 


7 TJ + 0 

ri c am 
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block 


(10.9) 


Sign conventions for the induction in Section 7 are opposite that in Figure 17, 
requiring the new notation 

a L = - and a T - -[W^ (10.10) 


so that the downwash vector is computed from 

W - -[0J - [$ 1 + [6 1 + [<f> 1 - [W T ] (10.11) 

ft 1 B 1 cam J ji L J ft lv block J /i 1 T J /i v ' 

A table in the computer printout lists the air angle accounting In these terms with 
one entry for each control point. 

There are other boundaries in the problem that should be mentioned at this 
point. The boundary at infinite radius has been taken care of in Section 3 by 
allowing only outgoing waves in the solution. Thus, we are guaranteed no reflec- 
tions or interference from any artificial outer boundary. Also in the wake as x-*» 
only convected vortex sheets appear with no possibility for reflected waves or 
influence back at the rotor. However, the boundary condition at the nacelle center 
body has been neglected for the current solution. Since the primary motivation for 
this work was to develop a noise prediction method, it was felt that errors in blade 
loading in the root area would not have an important influence on radiation predic- 
tions. Furthermore, blade root flow is usually dominated by transonic non-linear 
and viscous effects not covered by the present linear theory. 
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Induction Caused by Thickness 


In lifting surface theory for isolated wings, the flow field induced by the 
thickness effect can be ignored from symmetry considerations. For propellers, 
however, this symmetry no longer applies because of multiple blades and twist. 
Thickness at one radius can induce angle of attack at another radius. The pressure 
field for this effect was derived in Section 3 and can be converted into velocity 
potential and downwash following the procedures used in Section 7 for the loading 
effect. The derivation was given in Reference 15 in a slightly different notation 
and is so similar to the loading derivation of this report that it is not given 
here. Results can be written down almost by inspection. 

It turns out that the thickness effect is fairly weak, typically inducing flow 
angles of a degree or less. Thus, blade sections do not need to be represented very 
accurately. For expediency, we approximate all blade sections using a parabolic 
thickness distribution 


H(X) = 1 - (2X) Z 

and a thickness to chord ratio at each radius given by 

^ max thickness 
° chord 

The transform of this source is required by Equation 3.66 

w - sZ H(x) eiV “ 


where the wavenumber 


v = 2wa d 
O q 


The integral can be performed analytically with the result 
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The corresponding induced angle at the control points is 
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( 10 . 12 ) 

(10.13) 

(10.14) 

(10.15) 

(10.16) 
(10.17) 


This is in a form very much like the results of Equation 9.28, 9.29, and 9.36 for 
the loading induction and, regarding its evaluation, many of the remarks given in 
Section 9 apply. X arises from averaging over a small region near the control 
points and has the same meaning as in Equation 9.40 with q=0. It can be seen that 
the source transform for loading, given by X 1 in Equation 9.42, is replaced by V 

As a final remark, we point out that the parabolic thickness approximation used 
in this section for airload prediction is not used for noise calculations, where 
thickness may be the dominant effect. 
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SECTION 11 

AERODYNAMIC PERFORMANCE 


Propeller performance is measured in terms of thrust, power, and efficiency. 

In general, efficiency is power output divided by power input or 

_ THRUST * VELOCITY , n .. 

1 SHAFT POWER (11.1) 


Shaft power is an unambiguous quantity for this application but thrust is more 
complicated because of nacelle blockage effects. Apparent thrust relates to blade 
forces only whereas net thrust includes the penalty for nacelle forces associated 
with the propeller slipstream. Because this report deals only with blade forces, 
the efficiency computed herein is what the industry knows as "apparent efficiency" 
Thrust and power are used in coefficient form 


c T 


THRUST 
p N 2 D 4 

O 


Thrust coefficient 


and 


r POWER n „ . . 

C p — 3—5 , Power coefficient 

p N D 

o 

Combination of these 3 equations leads to 


> 7 



( 11 - 2 ) 


(11.3) 


(11.4) 


In previous sections treating blade forces, only lift has been considered. 

Lift, in the lifting surface context, is that force component perpendicular to the 
local direction of advance at any blade radius. If this were the only force acting, 
then the propeller would have unit efficiency because there would be no force compo- 
nent in the direction of motion. Drag is considered to have 3 components which can 
be discussed with reference to Figure 9. Wave drag is related to energy radiating 
through the surface S 2 and was treated in Section 6 via a sound power analysis. 
Vortex drag is related to the energy left in the wake and passing through S 3 . 

This can be computed via an integral over the vortex sheets. Finally, viscous drag 
at the blade surfaces is treated in strip fashion by the use of 2D airfoil tables. 

In the paragraphs below it is shown how thrust and power are related to lift and 
drag, how vortex drag is computed, how performance is computed on a linear basis, 
and then how this is corrected for non-linearity associated with finite axial induc- 
tion. 

Force Resolution 

fLe relation between thrust and torque components and lift and drag components 
can be deduced from the velocity triangles in Figure 18. The thrust force per unit 
radius is dT/Dr and the tangential force per unit radius is rdQ/dr where Q is 

torque. The corresponding lift and drag forces per unit radius are (1/2 )p U 2 bC. and 

q o L 

(l/2)pU 0 bC D . From consideration of similar triangles it can be found immediately 
that 
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(11.5) 
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r S - K U2b [ C L i + C D “] 


dr 2'o 

In terms of the coefficients defined above these become 

dC T 7T Z B<7B n 


dz 


4 a 2 


-(azC L - Cp) 


( 11 . 6 ) 


(11.7) 


and 


dC p 

dz 


7T 3 BctB d 

4a 3 


(C L + azC D ) 


( 11 . 8 ) 


No physics has yet entered the analysis; the above expressions simply represent a 
rotation of coordinates according to Figure 18. The components of C D are dis- 
cussed below. 


Vortex Drag 

Two dimensional airfoils in steady, subsonic, inviscid flow have no drag and 
trail no vortex systems. However, in 3D flow, lift falls off at the tips and a 
vortex system results. The drag associated with the energy of the vortices is known 
as vortex drag. Other names for the same phenomenon are induced drag and drag due 
to lift. Our analysis of vortex drag is an extension of the far wake method applied 
to wings in the book by Ashley and Landahl (ref. 27). This method is attractive 
because it applies equally to subsonic and supersonic blade section speeds and deals 
only with wake formulas, which are incompressible. In Figure 9, which was adapted 
from Reference 27, the flow through S 3 has a kinetic energy that, in principle, 
could be computed from the volume integral of the kinetic energy density. However, 
in the reference it is shown how to convert this volume integral to an integral over 
the surface of the vortex sheets 

K.E. = - P -f JJ A0 w dS ( U - 9 > 

where the jump in potential A ,<f> is equal to the circulation at corresponding span- 
wise location on the lifting surface and w is the induced downwash. The element 
of surface area in the notation of Section 7 is drdy where 7 is the surface 
coordinate in the advance direction in Figure 4 (not circulation) . The energy per 
unit area of the vortex sheet can then be written 


d[d(K.E. )/dr] _ . P _Q m » (11.10) 

d7 2 

To relate this to force on the local blade sections , note that an increment of 
kinetic energy can be related to an increment of work by the drag force as follows 

d(K.E.) - dF D dy (11. H) 
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or in per unit span terms 


d(K .E.) dF^ 

dr dr 07 


Thus , we have 


dF D p 


f- A<f> w 


dr 2 

This can be expressed in terms of the coefficient of vortex drag via 

dF P _ r 1 T,2 v 

dr C d 2 P o U b 


( 11 . 12 ) 


(11.13) 


(11.14) 


When these formulas are combined with the expression for section circulati 


ion 


, bUC. 

W 


the result is 


(I) --HS) 

VD 


far wake 


(11.15) 


(11.16) 


This is exactly the same expression that is well known in wing theory. It was 
rederived here to illustrate its origin and to demonstrate that it retains its form 
in the propeller context. 


An expression for the far wake downwash has already been derived in Section 7 
and will be adapted for vortex drag calculation. Equation 7.57 evaluated on the 
vortex sheets Is 



far wake 



S(z o ) 


l 


m=l 




(11.17) 


with n— mB . In the far wake, the source term is 

S(z °> ” ° 0 / AC P (7 0 .z 0 ) (11.18) 

which integrates to 


S ( Z o) “ 2 Vo C L (11.19) 

Combination of these expressions leads to the working equation for the coefficient 
of vortex drag 



-B 

2ttz 




B dC l 


1 I n (naz < )K sl (na z > )d z 0 

m=l 


( 11 . 20 ) 


The double bar on the integral sign signifies the finite part of an infinite inte- 
gral or the Mangier principal value (ref. 27). Methods to deal with this were 
discussed in Section 7. 
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Performance Calculation 

All of the necessary ingredients for computation of performance on a linear 
basis have now been developed. The procedure works as follows. The blade pressure 
distribution is calculated using the inversion method described in Section 10. 
Pressure distributions are integrated chordwise to find the spanwise distributions 
of lift coefficient. The C L 's are then used with the above formulas to find 
vortex drag. C L 's are also used with 2D airfoil tables to look up viscous (or 
profile) drag coefficients. This is done for the proper camber, thickness, and Mach 
number with corrections based on simple sweep theory. These 2 components of drag 
are used with Equations 11.7 and 11.8 to find the radial thrust and power distribu- 
tions which, in turn, are integrated to determine the overall thrust and power 
coefficients from 


dC T 

dz^ 


dz„ 


and 


dCp 

dz„ 


dz„ 


( 11 . 21 ) 


Then apparent efficiency in computed from Equation 11.4. If Mach number is high 
enough, the coefficient of sound power is computed using the formula in Section 6 
and added to the denominator in Equation 11.4 to account for the associated effi- 
ciency loss. 

A correction scheme has been developed to deal with the non-linearity caused by 
finite axial induction and is described in the remainder of this section. This 
aspect of non-linearity was discussed at the end of Section 8 and was shown to be 
related to the mass flux through the propeller disc. Evaluation of the circumferen- 
tially averaged flow showed that the relation between thrust and axial induced 
velocity did not satisfy the momentum equation. However, by adding the induction to 
the flight speed, momentum could be satisfied in the average sense. The problem was 
that the axial induction is not known a priori so that an iteration procedure must 
be used. The key to the method is to recall that in Section 7 the factor l/V was 
moved outside an integral, recognizing that a correction would have to be made 
later. V is the velocity that determines the mass flux and Is easy to correct 
because is appears in a simple factor outside the Integral representing the kernel 
function. 

The iteration procedure is illustrated schematically in Figure 19. It is based 
on satisfying the momentum relation given by Equation 8.32 which can be cast in 
terms of radial distribution of thrust coefficient as follows 


dC T 


v 0 " 7 r 3 z(l+v x /V 0 ) dz 


( 11 . 22 ) 


v y is the axial induced velocity at the propeller disc. The kernel function 
matrix is first computed at the actual propeller advance ratio so that the corre- 
sponding V is the flight speed V 0 . Then performance is computed using the 
scheme just described. Part of the performance calculation is to determine the 
thrust distribution dC T /dz. This is used to find the momentum related induction 
from Equation 11.22 as a function of radius but on a circumferentially averaged 
basis. The induction is used to correct the kernel matrix in blocks at each control 
point radius by multiplying by the factor V Q /(V o +v xd ) where v xd is the axial induction 
at the mid- load point of the disc. The corrected matrix is inverted again and the 
loop is repeated until convergence is achieved. 
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This procedure makes a great improvement in power prediction at takeoff speeds 
and powers. At very high loading, it overpredicts power and needs to be refined by 
recognizing the effect of increased dynamic pressure in producing blade forces and 
in including the loading associated with vortex flow in the momentum balance. 
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SECTION 12 

FORMULAS FOR WAKE PREDICTIONS 


This section presents formulas that have been programmed for the 3 components 
of the velocity disturbance caused by blade thickness and steady loading. Since the 
velocity potential was derived in Section 7, all that remains are some straightfor- 
ward manipulations to arrive at the working formulas, which are listed below. 

Volume III of this report presents sample calculations and comparisons with test 
data. 

In Section 7, we not only derived the velocity potential, Equation 7.26, but 
differentiated it in the axial and tangential directions to find the corresponding 
velocity components. These were used with Equation 7.28, which is repeated here to 
find the downwash component 


w 13$ ir^ 3x + M Ml (12 1) 

U “ U U L 3x d<t> d£ J 

normalized by the local section speed. The corresponding formula for the component 
parallel to the advance direction is 
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and the radial component 

is s 

imply 











V 

1 

d$ 



(12, 

.3) 




V 

" V 

d r 





where this is normalized by the flight speed. The partial derivatives for the chain 
ru le can be obtained directly from Equations 3.11a and 3.12a, after recalling the 
change in sign convention for x from positive in the flight direction in Section 3 
to positive downstream in this section. 

When the indicated manipulations have been performed, we arrive at the follow- 
ing formulas . 


m- + .«) 

(in - wr + + 1 **'- •-) 


(12.4) 


(12.5) 


( 12 . 6 ) 


Where ct sw is the steady wake term discussed at some length in Section 7 and (v/V) 5W is 
the corresponding wake term for the radial velocity. Note that there is no wake 
term in the u component. In applying these formulas, the program recognizes that 
the 2 wake terms exist only downstream of the rotor inside z— 1 . Elsewhere, they are 
equal to 0 , 
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The nearfield wake harmonic components are given by 
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( 12 . 8 ) 
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(12.9) 


The notation in these formulas has the same meaning as in Section 7 with the 
exception that X results from a Gaussian average over the field point in the 
axial direction rather than in the chordwise direction. It is given by 


X = exp 



(12.10) 


where A represents the width of the Gaussian function at its 1/2 amplitude points 
normalized by r T . 

The thickness and loading sources are represented by their transforms and 
as defined in Section 3. In the current version of the program, the section thick- 
ness distribution is approximated with a biconvex parabolic shape, whose transform 
was given in Section 10 as Equation 10.16. The loading transform is defined as 
follows 


- 1/2 

tf L (k x ) = I AC p (X)e ikxX dX (12.11) 

** ** 1/2 

This is similar to the definition used in the noise formulas in Section 3 except 
that Equation 3.74 was the transform of the loading shape only whereas Equation 
12.11 includes the lift coefficient. In the program, is obtained by operating 

on the loading vector L v that results from a performance calculation. 

The nearfield wake harmonics given by Equations 12.7 to 12.9 are non- singular 
upstream of the blade leading edges and downstream of the trailing edges. There- 
fore, the radial integrals do not have to be divided into singular and non -singular 
ranges as they did for kernel function calculations. 
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SECTION 13 

STUDIES WITH 3D UNSTEADY LIFT RESPONSE THEORY 


In Sections 3 and 7 above, the unsteady lifting surface integral equation was 
derived relating fluctuating downwash at the blades to unsteady lift pressure. 
Section 9 explained the system for discretizing the equation and Section 10 showed 
how to invert it. This capability has been coded and delivered to NASA. Input to 
the program is described in Volume IV of this report and some applications to real 
geometry are given in Volume III. In this section, we show the results of some 
cases that were run to verify correct behavior of the theory by comparison with 
well-known 2D theories and then to explore how 2D and 3D results differ. The cases 
run are intentionally academic (constant chord blades) so as to concentrate on 3 
dimensional and compressibility effects. 

The following paragraphs describe the program's capability and limitations and 
interpretation of the input and output. Then cases for response to gusts and blade 
plunging motion in incompressible flow are compared with the Sears and Theodorsen 2D 
theories. Finally, a 10 per revolution frequency case is studied for effects of 
compressibility, sweep, and comparison with Amiet's 2D compressible gust response 
theory (ref. 34). 


Program Cap ability and Limitations 

As currently programmed, the code can deal rigorously with the following varia 
tions in geometry and operating conditions. 

- The program predicts the unsteady lift pressure distribution caused by 
either non-uniform inflow, blade vibration, or any combination of the 
two. Any harmonic motion involving twist, camber and bending oscillation 
can be treated, 

- The program deals with one frequency at a time. The user specifies the 
frequency in terms of P order (integer or fractional). Any physically 
permissible interblade phase angle can be run (input is in terms of 
number of nodal diameters) . 

- The theory is 3 dimensional with unshrouded blade tips and deals rigo- 
rously with blade interference, sweep, taper, and Mach number variations 
along the span (including supersonic tips). 

- Bow and trailing shocks and their interaction with adjacent blades are 
represented as in linear supersonic wing theory. 

There are some limitations to the code in its present configuration and also some 
limitations to the basic theory that should be noted. 

- The theory is strictly linear. This means that the unsteady distur- 
bances are not coupled to the steady induction except possibly in the 
circumferentially averaged sense by running at an adjusted advance ratio. 
Moment oscillations caused by a mid-chord oscillating shock would not be 
represented. 
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- The unsteady portion of the program does not deal with vortex loading. 
This may be a significant shortcoming because the small extra load due to 
a tip edge vortex could cause significant bending moments at the root. 

Frequencies are limited to about 20P (20 times the propeller rotation 
frequency) for the cases examined so far. This is not an intrinsic 
limitation and can probably be extended without much effort. 


General Input /Output 

Input falls into the categories of blade geometry, mean operating conditions, 
calculation parameters, and unsteady boundary conditions. This last category is 
treated in more detail under the subsection on Details of Panels and Input/Output 
and in the sections showing sample calculations. The major input items are 

- Geometry ^Number of blades 

*Radial distributions of Mid- Chord 

Alignment, Face Alignment 
Diameter, Hub/tip ratio, and 
Chord/diameter 

- Operating Conditions: ^Advance ratio 

^Flight Mach number 

Calculation Parameters: ^Number and location of control points 

Boundary Conditions: *P order (frequency/prop rotation 

freq), q 

^Number of nodal diameters 
^Unsteady downwash velocity vectors 
(complex) at each control point 

The program prints out the complex pressure coefficient AC p , the centers of 
pressure for the real and imaginary parts of AC p , and lift coefficient C L . To 
recover a dimensional lift pressure, note that 

AP - I p o U z AC p (13.1) 

where is the ambient density and U is the section relative velocity such that 

U 2 - V 2 + z z V 2 (13.2) 

V x anc * are *-he flight speed and tip rotational speed and z is radius ratio. 

The actual program representation is 

U 2 = M^c z [ 1 + a z z z ] (13.3) 

where M* is the flight Mach number, c Q is the ambient speed of sound, and a - 
n/J . Similarly, lift per unit span is 

L = ^ 0 uZbC L (13.4) 

where b is section chord. 
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Details on Panels and Input/Output 


Figure 14 shows the load paneling system. The blade is divided into spanwise 
panels which are tapered at constant percent chord. The number of panels is speci- 
fied by the variable NCP. The control points, where the flow tangency boundary 
condition is enforced, are at the midpoint of the panels in chordwise arrays at 
radii 7 . i . The index system for the control points is the same as in reference 2, 
namely that i=l,2, ... NSM counts the radial locations and m = 1,2, ... NCP 
counts the chordwise positions. The current limits on NSM and NCP are 10 and 20, 
respectively. The downwash velocities are input as the vector where 

H = (i-l)NCP + m (13.5) 


and the velocities are normalized by the local section speed U. 

The loading elements have spanwise and chordwise variations as follows. In the 
spanwise direction the variation is given by a series of shape functions or spanwise 
modes as shown at the top right in Figure 14. For the chordwise variation, a 
series of overlapping triangular elements is used with some special treatment at the 
leading and trailing edges as shown at the lower right in Figure 14. At the leading 
edge, the user can choose the triangular element shown with the dashes (NBAR0PT=1) 
or the prefered singular element (NBAR0PT=2) shown by the solid line. The latter 
element is rigged so that, with the adjacent triangular element, it gives the cor- 
rect leading edge behavior for 2D cases. At the trailing edge the basic element 
goes to 0 for subsonic trailing edges [M rel cos (trailing edge sweep angle) <0] satis- 
fying the Kutta condition. For supersonic edges a finite jump is allowed since the 
Kutta condition doesn't apply. There is a blending at the sonic radius between the 
2 types of element. 

The loading elements are simply shape functions. To find the actual loading 
they are multiplied by coefficients which must be determined by matrix inversion 
so as to satisfy the boundary conditions. The index u is defined by 

i/ - (j-l)NCP + h (13.6) 

The associated loading is given by 

<13 ' 7> 

which is continuous in the chordwise direction and continuous and smooth in the 
spanwise direction. AC p is printed and plotted with the computer run and is sent 
as a table to the noise program for unsteady loading input. 

The loading AC p and disturbance vector W are complex quantities. To 
interpret them as physical quantities, note that the downwash at any control point 
can be written 


where q is the P order 


g - W e iqnt = (w r + iw A ) e" iqflt (13.8) 

(i.e. frequency normalized by the shaft rotation frequency). 
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The physical quantity is the real part of this, namely, 

w r cos (qflt) + sin(qOt) (13.9) 

Similarly, the physical pressure coefficient and lift coefficient are the real parts 
of AC p and C L . For example, for the pressure the real part is 

AC p cos(qfit) + AC P sin(qflt) (13.10) 

real imag 


Verification of Correct Response to Gusts and Plunging Motion 

The first step in verification of the new blade pressure response theory is to 
compare predictions with results from well-known 2D lift response functions. This 
serves to check for any coding errors as well as to establish how the 3D flow dif- 
fers from 2D strip behavior. The classic 2D unsteady lift response theories are 
Sears', in which the airfoil is rigid and the inflow contains a sinusoidal upwash 
gust field, and Theodorsen' s , in which the inflow is steady and the airfoil executes 
a sinsoidal plunging motion. Both theories apply to linear, incompressible flow 
about zero thickness flat plate airfoils with no mean loading coupling. The text by 
Bisplinghof f , Ashley, and Halfman (Ref. 3) provides good background on the deriva- 
tion and use of these theories. In each of the following 2 subsections, the results 
of the 2D theories are summarized first and then compared with 3D predictions. 


Comparison with Sears Lift Response Function - Sears' results are outlined in Figure 
20. The airfoil with chord b is fixed at the origin in a flow field with mean 
velocity U. Superimposed on U, and convected with it, is the upwash gust w 
described by 


w = w o e iu(t " X/U) 


(13.11) 


where the frequency f is w/27r. Here we use Sears' original convention with a + 
sign for time in the exponential (+iwt) . Because x = Ut is a constant phase point, 
it is clear the equation 13.11 represents a gust convected with the free stream. 

The standard reduced frequency parameter is based on semi -chord 


wb/2 

V 


(13.12) 


Since wavelength A = U/f and w -= 27rf, the reduced frequency can also be inter- 
preted in terms of chord to wavelength ratio: 

k 0 - (i3.i3) 

For a small steady upwash w Q , the incidence angle is a = w Q /U and the lift coeffi- 
cient is C L - 2na. Sears' result for the unsteady case with a = w/U is 

C L = 2?raS(k o ) (13.14) 
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where S(k Q ) is the Sears function which describes the deviation in amplitude and 
phase from the steady theory. This is plotted in Figure 20 with reduced frequency 
as the parameter. It can be seen that the function goes to 1 at 0 frequency, as 
required, and diminishes in amplitude at higher frequencies. The phase first leads 
the input then lags as frequency is increased. 

Input for the gust cases was prepared as follows. With the - iwt sign conven- 
tion used for the 3D analysis, a sinusoial gust is described by 

w = w o e i,jC7/lJ ‘ (13.15) 


where 7 is a streamwise coordinate measured in the helicoidal advance surface at 
constant radius. 7 corresponds to x in the 2D discussion above. The input vector 
W Is obtained by dropping the exp(-iwt) and normalizing by the local section speed 
U to get the upwash angle : 



^ iu-7/U 

U 


(13.16) 


w o /U was set to 1 for a 1° gust amplitude. We substitute w - qft, x = bX- , b = 2r T B D , 
flr T /V - a, and U/V * o with the result 


W = cos|^ 


f 2 q aB D 




, f 2 q aB D 


sin 


*5] 


(13.17) 


which is to be evaluated at the non-dimensional chordwise locations 


X- - -0.5 + (n - 0 . 5)/NCP (13.18) 

where h runs from 1 to NCP. 10 control points were arrayed across the 0.7 radius 
line so that NCP = 10 and NSM « 1. Cases for the geometry shown in Figure 21 were 
calculated for P orders from 0 to 20 with the results shown in comparison with the 
Sears function. To interpret the P order in terms of reduced frequency, the same 
substitutions used for Equation 13.17 can be used with the definition of k Q with 
the result 


k 


o 


q aB D 

O 


(13.19) 


For plotting on Figure 21, the lift output from the 3D lifting surface program was 
normalized by 27ra, which amounts to dividing by 0.1097 after converting to 
radians. Furthermore, the complex conjugate is plotted because of the difference in 
conventions on the sign of iwt in the exponentials of Equations 13.8 and 13.11. 

The comparison of 2D and 3D results in Figure 21 can be interpreted as follows. 
In the limit as frequency goes to 0, the steady result for a twisted flat plate 
blade is reached. The lift is about 1/2 of the 2D result because of the 3D induc- 
tion. This lift reduction corresponds to the aspect ratio effect, well known in 
steady wing theory. At the high end of the frequency range, 20P, the 3D lift 
approaches the 2D value (at least away from the ends of the blade) . This effect is 
also well known in unsteady wing theory (Ref. 35) and is a very satisfying check of 
the 3D theory. 
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Since the high and low frequency ends of the 3D curve in Figure 21 could have 
been predicted a priori, the valuable part of the result is in the detail of the 
transition from steady 3D behavior at low frequencies to unsteady 2D behavior at 
high frequencies. The surprise was that at IP there is almost no effect of 
unsteadiness. Thus, for this incompressible, unswept case, the 3D lifting surface 
theory theory gives results which should depart only minutely from a quasi -steady 
lifting line calculation. This is an important discovery and answers questions 
regarding unsteady lift response effects in IP calculations that have been raised 
over the past 30 years. In the past at Hamilton Standard, various analysts have 
attempted to account for unsteadiness by doing a multi-azimuth steady lifting line 
calculation and then adjusting the amplitude and phase in a strip-wise sense by 
using the Sears function. Figure 21 shows clearly that this is not correct and 
that, fortuitously, the quasi-steady calculations should be adequate for IP and 2P, 
at least for this unswept case. 

This conclusion regarding unsteady effects can at this point only be applied 
rigorously to the cases run so far and shown in Figure 21. However, based on other 
experience, we would expect it to apply also up to high subsonic tip speeds for 
straight blades. It is recommended, however, that this be explored further by 
trying more cases so that rules of thumb can be developed on the limits where quasi 
steady multi -azimuth calculations can considered accurate. 

The final portion of the Sears function check is a comparison of 2D results 
with 3D unsteady pressure distributions calculated by the program. Sears' theory 
modifies the flat plate steady results in amplitude and phase only. The shape of 
the pressure distribution stays the same: 


ac p “ 4 “\F5? S(k Q ) (13.20) 

Comparisons are given in Figure 22. As might be expected, the 2D and 3D pressure 
agree closely at the high frequencies (20P) where the lift coefficients matched. 

This applies to both the real and imaginary parts (in phase and out of phase) and to 
the detailed behavior near the leading edge. At lower frequencies the pressures 
depart from the 2D results in a way that is consistent with the lift predictions. 

The distributions at the upper right in Figure 22 have the shape given by the 
radical in Equation 13.20 with the characteristic square root singularity at the 
leading edge and the zero (Kutta condition) at the trailing edge. Later in this 
section, distributions with these characteristics are referred to as "flat plate 
like . " 

Comparison with Theodorsen Theory for Plunging Airfoils - Results of Theodor - 
sen's theory for the plunging airfoil problem are outlined in Figure 23, which 
includes several parameters defined in the Sears discussion above. Here the 
upstream airflow is steady and the airfoil is executing an oscillating vertical 
motion given by 


h - h 0 e ivt 

The effective unsteady angle of attack caused by this motion is 


(13.21) 


a = h/U 


(13.22) 
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Theodorsen's theory deals with various kinds of motion, but his result for unsteady 
lift in the plunging problem is" 

C L = 2na [ C (k Q ) + ik c /2] (13.23) 

where C(k Q ) is the Theodorsen function as given in Ref. 3. C L /2na is plotted in 
Figure 23 with reduced frequency k 0 as the parameter. At the low frequency end of 
the curve, the value is 1, corresponding to steady, 2D airfoil lift. For increasing 
frequency, the lift first leads the motion and then lags. At high frequency the 
lift continues to increase because of the inertial mass effect. 

For verification of the 3D theory, the same propeller and same mean conditions 
were used as with the Sears theory above. Since the motion is constant along the 
chord, the input unsteady boundary condition for the program is simply 1 ' s at all 10 
chordwise control points corresponding to a 1° effective angle of attack due to 
the airfoil motion. Again, for comparison with Figure 24, the output lift coeffi- 
cient from the program is divided by 0.1093 to normalize by 2na and the complex 
conjugate is taken. As with the Sears theory, the 3D results exhibit an aspect 
ratio effect at low frequency and approach the 2D curve at high frequency. This 
general behavior is to be expected from wing experience (Ref. 35). 

Theodorsen' s unsteady pressure distribution for the plunging problem can be 
deduced from Reference 3 : 

AC p - 4a C(k 0 ) + i 2k 0 <Jx(l-x) j (13.24) 

Figure 25 shows the comparison between 2D unsteady pressure predictions from Equa- 
tion 13.24 and predictions from the 3D unsteady theory. As with the Sears' theory, 
3D results approach the 2D unsteady theory at high frequency and differ signifi- 
cantly at low frequency. 


10P Cases for Noise Analysis at Takeoff Conditions 


The test cases presented above for checkout against 2D incompressible theory 
were run with control points at only one radial station. We now examine a series of 
3 cases run with 10 control points across the chord at each of 7 radii to find how 
the 3D flow compares with the 2D strip theory currently in use for noise predic- 
tions. These cases include compressibility at Mach numbers and advance ratios 
typical of takeoff. To minimize the number of effects in the study, the blade again 
has constant chord as shown in Figure 26. There are 5 blades and the interblade 
phase angle is 0 as would be the case for a 5x5 CRP. The frequency is 10P, which is 
the fundamental interaction frequency for counter-rotation interaction. Flight Mach 
number is 0.26 and the gust amplitude is 1° from root to tip in all cases. 

The 10P series is shown in Figures 26 to 29. Figure 26 is a base case with no 
sweep and a J for a low tip relative speed (0.64 Mach). The shape of the pressure 
distribution follows the typical flat plate behavior expected for low Mach number 
and frequency. Amplitude and phase effects will be examined shortly. 

The next case considered is for swept blades. To evaluate the boundary condi- 
tions for program input, we consider the same gust given by Equation 13.15. With 
the blade swept back in the gust field by an amount MCA, the blade sees the gust at 
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a time lagged by At — MCA/U. Thus, t in Equation 13.15 can be replaced by t-At. 
The associated phase lag is t*>At, or in non-dimensional terms 


A (f> 


2q a MCA 
a D 


which changes the input vector from Equation 13.17 to 


(13.25) 


W = w 0 



+ 


2qa MCA 1 

° D J 


+ i sin 


[ 


2qaB D 

a 


X- 


n 


+ 


2qa mca "I 1 

° D J J 


(13.26) 


Figure 27 shows results for a uniform 45° sweep. There are rapid phase variations 
along the span; however, both the in phase and out of phase distributions have flat 
plate like behavior. At radii where the pressure is in transition from positive to 
negative, the unusual shapes are probably simply an indication that the pressure 
doesn't change sign simultaneously along the chord. 

Results from the previous 2 figures are summarized in terms of lift amplitude 
and phase in Figure 28. At the top, the reduced frequency and section relative Mach 
number distributions with radius are shown for reference. The middle panel shows 
the magnitude of the lift coefficient compared with Amiet's compressible lift 
response function (Ref. 34) used in a strip sense. It can be seen that the strip 
method is a good approximation for the straight blade over all of the span except 
for the tip where some 3D relief is found. The results for the swept case show some 
lift reduction in the mid-blade area, which, of course, is a benefit for noise. The 
phase comparison is given at the bottom of Figure 28. For the swept blade compari- 
son, Amiet's phase is lagged by an amount corresponding to how far the blade sec- 
tions are moved back in the gust field. For both the straight and swept cases, 
Amiet's theory used in strip fashion agrees well with the 3D theory in the mid-blade 
area but deteriorates toward the tip. 


Figure 29 shows results for a straight blade at higher tip Mach number. Here 
the pressure distributions depart significantly from the simple Sears shape. At the 
tip, there are some oscillations along the chord, which also occur in the 2D com- 
pressible theory although no quantitative comparisons have yet been made. In the 
root area there are also some unusual shapes, presumably caused by cascade effects. 

The conclusion from these 3 cases is that the Hamilton Standard strip methods 
could be adequate at the 10P frequency for low tip Mach number. Higher frequencies 
and Mach numbers need further exploration. 
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SECTION 14 

CONCLUDING REMARKS 

Theoretical equations have been developed for computation of noise and aerody- 
namics of propellers associated with blade loading and thickness. Derivations were 
based on the acoustic analogy and the acceleration potential method and therefore 
are basically linear. However, in the case of propeller steady performance, non- 
linearity associated with the axial mass flux is crucial and a correction method was 
presented. This essentially linearizes the induced velocity about a mean flow 
through the rotor that satisfies the axial momentum equation on a circumferentially 
averaged basis. 

To a large extent this volume extends and formalizes theory for aerodynamics 
and noise that has been presented by the author in various technical papers. Wher- 
ever possible, the relationship between formulas developed herein and previously 
published work is shown. In some cases, the formulas are equivalent but with a 
different notation; in other cases, the formulas presented herein are more general 
than those in previously published work but reduce to the earlier theories when the 
correct restrictions are applied. 

This volume is mainly theoretical and presents little in the way of computed 
results. However, Section 13 does show some calculations verifying correct behavior 
of the 3D unsteady lift theory in comparison with traditional 2D methods. Verifica- 
tion of the theory by comparison with test data is the subject of Volume III of this 
report. Volume IV is a user's manual for the computer code. 
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SECTION 16 

LIST OF SYMBOLS 


There are a few symbols not included in this list that are used only within a page 
or so of the point where they are defined. 

a — tip rotational speed/flight speed — tt/J 

b - blade section chord measured at constant radius - see Figure 3.1 
c Q — ambient speed of sound 
f = body force in wave equation 

f - j th component of force/unit volume on fluid 
f « j th component of force/unit area on blade 

h = section thickness distribution in Section 3, camber distribution in Section 
10; section vertical displacement in Section 13 

x — 4-~T, also index for control points in radial direction when used as a 
subscript 

j - index for radial mode shapes , 

k - number of nodal diameters in unsteady loading pattern - related to interblade 
phase angle - see discussion with Equation 3.80 

k Q = reduced frequency - Equation 13.12ff 

k x - normalized axial wavenumber - Equation 3.102 

ky “ wavenumber for sound calculation - defined in various contexts in Sections 
4 and 5 

k 7 = chordwise wavenumber - Equation 3 . 9 

m - harmonic index, nB for steady loading, nB-k for unsteady loading 

p ** disturbance pressure 

p^ - pressure field of rear rotor 

q - source (monopole) strength in Section 3, unsteady loading order elsewhere. 
Loading frequency is qQ/2n 

r — radius in cylindrical coordinates, radiation distance in Sections 5 and 6 
only 

r T - tip radius 

t - time 
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t b — ratio of section max thickness to chord 


u = component of disturbance velocity parallel to the local section advance 
direction 

v = component of disturbance velocity in radial direction 

w — component of disturbance velocity normal to advance direction - the downwash 

(negative down) 

x « distance forward of pitch change axis in acoustic sections (Sections 3-6), 
distance downstream of pitch change axis in aerodynamic sections (Sections 7 
- 12 ) 

“ x coordinate for evaluating kernel function elements 
Ax * axial distance between load panel and field point 

y = spanwise distance in wing discussions, sideline distance to field point in 
noise formulas 

z = radius ratio for field point or observer position 

z Q — radius ratio for source integration 

B — number of blades 

B d - b/D, chord to diameter ratio 

C L — section lift coefficient 

C D * section drag coefficient 

Cj^ — coefficient of radial tip force - Equation 5.44ff 
C- = load element shape functions - see Figure 9.1 
AC p * coefficient of lift pressure 
D = propeller diameter 

E — various exponentials - defined locally 

F — body force in momentum equation 

FA = face alignment measured at constant radius - see Figure 3.1 

G - Gaussian averaging function for control points - Equation 9.21 

H - Heavyside or unit step function - 0 for argument < 0, 1 for argument > 0. 

Also, normalized section chordwise thickness distribution - see Figure 3.3 

H^ 1} = J n + iY n , Hankel function of first kind 

— J n - iY n , Hankel function of second kind 
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I - modified Bessel function of first kind. 

n 

j - advance ratio, V/[ (revs/unit time)*D] 

J n - ordinary Bessel function of first kind 

^(w) — Bessel function product defined in Equation 7.27 

K - normalized radial wavenumber - Equation 3.55 

- modified Bessel function of second kind - Equation 3.17 

K — radial wavenumber - Equation 3.36 
r 


- ak 7 , axial wavenumber 

K^.etc. ^wavenumbers defined in Equations 7.17 - 7.20 
K — kernel function elements 

\IV 

L “ loading coefficient - see Equation 9.25 

- flight Mach number 

«= section helical Mach number 

M -= nr t /c 0 , tip rotational Mach number 

MCA = mid-chord alignment or sweep of section - see Figure 3.1 

NCP = number of loading panels 

NSM = number of spanwise mode shape functions 


Q 

R 

R o 

S 

T 

U 

V 

V 

X 


- torque in performance section 

«= source to observer distance in Green's function 

- see Equation 3.34 

*= source strength - see Equation 3.2. 

- Thrust 

_ blade section relative speed, - aV 

- flight speed 

- mean flow speed through propeller disc - a function of radius ratio 

- vector of downwash angles at control points 

- normalized chordwise distance, zero at mid-chord 
= ordinary Bessel function of second kind 
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induced flow angle - relative to section advance direction 


- Jl - K 

= coordinate at constant radius in section advance direction-see Figure 3.1 

- Dirac delta function (impulse function) 

- Kroneker delta, 1 for equal subscripts, 0 for unequal 
- angle in cylindrical coordinates to field point 

“ angle in cylindrical coordinates to source point in translating, but not 
rotating, system 

- angle to field point in blade -fixed system 

= phase angle associated with face alignment - Equation 3.104 
= phase angle associated with sweep - Equation 3.62 

= phase angle associated with axial position of field point - Equation 3.78 

“ *FA + * s 
“ ^FA + + 

“ see Equation 7.31 

* index for control points - see Equation 9.17 

* index for load elements - see Equation 9 . 24 

“ angle from flight direction for acoustic equations 
-* ambient density 

- apparent efficiency in propeller performance 

- 'Jl+a z z 2 , ratio of local blade section speed to flight speed at 
field point radius 

** Ji+ a2 z o 2 > ratio of local blade section speed to flight speed 
at source point radius 

= wavenumber integration variable for axial direction 

- (unsteady loading frequency) times 2n 

- distance normal to local section advance direction measured at constant 
radius - see Figure 3.1 

— source transforms defined in Section 3, with various subscripts 
- source time integration variable in Green's function 



$ = velocity potential 

r - 1 / r t 

r a - y o ^t 

II - sound power in Section 6 

n = propeller angular speed, 2n x shaft rotation frequency 
* - normalized versions of 

- loading source transform - Equations 3.72 and 3.69 

- combined loading and thickness transform - Equation 3.69 

UNDERBAR 

denotes vector quantity as in F 
OVERBAR 

denotes average quantity 
SUBSCRIPTS 

i - index for radial position of control points 
k - k th loading harmonic 
n - n th sound harmonic 

o - denotes source coordinate or properties 
T * denotes tip position 
h - denotes hub position 
< - smaller of r,r 0 or z,z o 

> - larger of r,r D or z,z o 

INTEGRALS 

£ -= denotes Mangier principal value for integral 
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Figure 1. Prop -Fan model SR3. 
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Figure 2 . 


Vortex flow visualization on a simulated marine propeller in a water 
tunnel. (from Reference 18) 
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Figure 3. Flow visualization results for CRPX1 at take off condition. Change in 
streak direction indicates vortex reattachment line. 
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Figure 6. Normalized section loading 
F x (X) and H(X) . 
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Figure 8 . Retarded and visual coordinates . 
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Figure 10. 


Passage of fluid particle past fixed wing and through rotating rotor. 
Particle misses airfoil in both cases. Particle path not deflected in 


linear theory. 
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11. Actuator disc induction velocity components. 
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Figure 12. Sample calculation of axial induced velocity for typical Prop-Fan geom- 
etry. Note developing jet and forward induction outside the tip radius 
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Figure 13. Axial induced velocity distribution for rotor from previous figure at 
radius ratio z *= 0.78. Flight Mach number « 0.72. 
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CONTINUOUS LOAD DISTRIBUTION 
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Figure 14. Load paneling system. 
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Figure 16. The 3 types of load panel. 
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Figure 17. Angles 


used to define 


boundary conditions . 
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Figure 18 . Resolution of forces for calculation of performance . 
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Figure 20. Airfoil/gust interaction 


problem for Sears Theory. 2D, Incompressible 




Figure 21. Propeller blade/gust interaction problem for 3D analysis. Incompressible 
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Figure 22. Unsteady pressure coefficients. 3D versus 2D. 
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Figure 25. Unsteady pressure coefficients. 2D versus 3D. 
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Figure 29. Gust interaction at 10P. Higher Mach number. 
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